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FOREWORD 


The growth of aeronautical engineering and science, 
particularly during recent years, has been more rapid 
than the means for disseminating the experimental and 
analytical knowledge acquired during that growth, or of 
presenting the results of much of the original work done 
by the Industry, Universities, and Research Establishments. 

Progress in aeronautics is dependent on the results 
of these researches and on the development of new and 


improved methods of analysis or experimental techniques, 


and the publication of such original work will not only 
stimulate and encourage the workers engaged in these 
activities, but should ensure that their work is more 
widely known and applied. 

Under the terms of its Charter, the Royal Aero- 
nautical Society ‘s charged with the duty of advancing 
Aeronautical Art, Science and Engineering and, among 
other things 

‘“‘to facilitate the exchange of information and ideas 
amongst the members of the Society and others.” 

In carrying out this task, it is the primary object 
of the Editorial Board of The Aeronautical Quarterly to 
encourage workers to submit papers describing new and 
original work, or papers reviewing the progress in some 
specialised field of activity and, through the medium of 
this new publication, to make the results of these 
researches available to all workers or design groups 
concerned with aviation. 


THE AERONAUTICAL QUARTERLY 


To ensure that the high standard set by the Council 
of the Society shall be reached and maintained, the 
Editorial Board has been fortunate in having the willing 
co-operation of the leading authorities in various branches 
of aeronautical engineering and science to act as referees 
on all papers submitted for publication. 


The Council of the Society was greatly encouraged 
to sponsor this new venture by the response it received 
from members of the Society to its enquiry as to the 
need of such a publication. The future of The Aeronautical 
Quarterly is dependent on the support it receives; it is 
hoped, therefore, that this support will be forthcoming 
from all sections of the Society, so that the success 
of the new venture is assured. 
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CONTROL REVERSAL EFFECTS ON 
SWEPT-BACK WINGS 


by 


HAYDN TEMPLETON, B.Sc., A.F.R.Ae.S.* 


(Structures Department, Royal Aircraft Establishment) 


SUMMARY 


Aileron reversal effects on swept-back wings in general and elevon reversal 
effects on tailless swept-back wings in particular are discussed on a non-mathematical 
basis, attention being confined to the orthodox flap type of control. The main 
purpose of the paper is to convey in the simplest terms possible a clear physical 
picture of the conditions producing loss of control power, emphasis being naturally 
laid upon the part played by structural wing distortion. Certain qualitative features 
relating to the two phenomena are also discussed. As a general introduction to the 
discussion on aileron reversal effects, the definition of “aileron power” in relation 
to the actual dynamic condition of rolling is described at some length. For elevon 
reversal effects on tailless aircraft the effect of wing flexibility on both “elevon 
power” and on trim in steady symmetric flight is considered. 


With the descriptive treatment adopted the analysis is of necessity broad and 
general but is designed to appeal to those not too familiar with the subject. The 
results of certain calculations on a hypothetical wing, which may be of interest, 
are included. A mathematical analysis for the quantitative estimation of both 
aileron and elevon reversal effects is given in the Appendix. 


1. WING DISTORTION 


Implicit in the study of control reversal effects is the consideration of how a 
wing distorts under the action of applied loads. The phenomenon of control 
reversal effects may in general terms be described as follows. A particular control 
is deflected in order to produce a certain aerodynamic effect, but the aerodynamic 


Paper received November 1947. 


* This paper is based on work done while the author was with Sir W. G. Armstrong-Whitworth 
Aircraft Ltd. 
[The Aeronautica! Quarterly, Vol. 1, May 1949) 
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loads resulting from the control deflection distort the flexible wing and the aero- 
dynamic effect of the wing distortion modifies the effect that would have been 
obtained from the control on a rigid wing. The loads which produce the wing 
distortion are finally the net result of the loads due to the control and the loads 
due to the wing distortion itself. In the calculation of aeroelastic effects two 
problems therefore present themselves: the determination of the elastic deformation 
of the wing structure under an applied loading system, and the determination of 
the aerodynamic loads for a given wing configuration. Both these are major 
problems in themselves, one structural and the other aerodynamic, but when 
combined as they are in aeroelastic phenomena, where wing distortion and aero- 
dynamic loads mutually affect each other, then some simplification is inevitable if 
a practicable solution is to be obtained. 


Methods used so far fall into two main categories. Either the semi-rigid 
principle is employed in which arbitrary modes of distortion are used and the 
solution based on the determination of the absolute displacements at some reference 
section; or the mode is correctly determined from the elastic characteristics, in 
which case an iterative process must be used. A combination of these two methods 
is described later in this paper. Whatever method is used, however, it seems that 
for a reasonably workable solution simplifying assumptions must be applied to 
bot the structural and aerodynamic sides of the problem. 


On the structural side it is a question of how the wing distorts under typical 
loading actions and of how the distortion should be represented in terms of its 
aerodynamic significance. As a preliminary it is convenient to establish certain 
useful definitions with regard to direction. The chord is measured parallel to the 
longitudinal axis of the aeroplane and a chordwise direction is therefore one 
parallel to this axis. Normal direction is that at right angles to the wing, or more 
specifically to the wing flexural axis (see Figs. 2 and 3). These are the definitions 
adopted here and are not of course necessarily universal. 


If the aerodynamic loads are expressed in terms of lift and pitching moment 
on the chordwise section, then it follows that any distortion which involves a change 
in incidence or shape of the chordwise section is significant aerodynamically. A 
plain change in incidence produces a change in lift at the aerodynamic centre: a 
change in shape is more complicated but can be represented as an effective change 
of camber, producing changes in both lift and moment at the aerodynamic centre. 
It should, incidentally, be stated here that as discussion is confined to steady or 
effectively steady states, only the displacements due to distortion, and not velocities 
or accelerations, are to be considered. Pure translation of the chordwise section in 
a direction normal to the air stream, for instance, has in the steady state no aero- 
dynamic effect, whereas velocity of translation in this direction has of course a direct 
aerodynamic effect resulting from the effective change of incidence due to the 
superimposed velocity. 
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WL: 
PAR 
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Fig. 1. 
Unswept wing with chordwise ribs. 


Considering first of all the case of an unswept wing with the usual arrangement 
of ribs running chordwise from the spar (Fig. 1), it is evident that as far as the 
spar itself is concerned only torsional deformation is significant, resulting in rotation 
of the chordwise section. Flexural deformation produces only a translation of the 
chordwise section and has therefore no effect. Deformation of the ribs may, in 
addition, produce both rotation and distortion (change of shape) of the chordwise 
section, and deformation of the skin between the ribs could have a similar effect. 
To take account of all these possible sources of distortion would be a formidable 
task and in the estimation of aileron reversal effects on unswept wings it has been 
the practice to consider only the torsional deformation of the spar, excluding the 
more local deformations of ribs and skin by what amounts in effect to an assumption 
of rigid closely-spaced ribs. 


Coming now to the swept-back wing, there are two different arrangements 
to be considered. In the first (Fig. 2) the ribs again run chordwise. If the assumption 
of rigid closely-spaced ribs is maintained there can be no distortion of the chordwise 
section, but as a purely geometrical consequence of the sweepback both torsional 
and flexural deformations of the spar produce rotation of the chordwise section. 
Relative to the centre say, torsional deformation results directly in a rotation # of 
the normal section which produces a rotation @ cosy of the chordwise section, 
y being the angle of sweepback. Flexural deformation, on the other hand, results 
directly in a change of slope ¥ along the flexural axis which produces a rotation 
u sin y of the chordwise section. Increasing the sweépback therefore respectively 
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Fig. 2. 


Swept wing with chordwise ribs. 


Fig. 3. 
Swept wing with normal ribs. 
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reduces and increases the contributions due to torsional and flexural deformation 
of the spar. These rotations of the chordwise section relative to the centre constitute 
an aerodynamic wash-in or wash-out depending upon the sign of the basic deforma- 
tion. They may be summarised as follows : — 


Basic deformation 


(relative to centre) Aerodynamic effect 


1. Nose-down twisting 
Wash-out 
2. Upward flexing 
Nose-up twistin 
Wash-in 
4. Downward flexing 


In the second arrangement to be considered (Fig. 3) the ribs run in the normal 
direction and with the same assumption of rigid ribs there is therefore no distortion 
of the normal section, from which it must be concluded that there is distortion of 
the chordwise section. This distortion is probably not important in the case of 
torsional deformation of the spar, but with flexural deformation it appears as a 
direct change of camber. The effect is illustrated in Fig. 4, which really does no 


CAMBER IN CHOROWISE 


‘UPWARD FLEXED 
DIRECTION. 


WING. 


Fig. 4. 
Camber change on swept wing. 
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more than point the well-known fact that in a curved surface there can be only one 
direction in which there is no curvature. The amount of the camber change is a 
function of the angle of sweepback and of the curvature of the flexural axis in 
the mode of flexural deformation: its sign is directly related to that of the flexural 
mode, as follows : — 


Upward flexing — negative camber change 


Downward flexing — positive camber change 


In addition to distortion, the chordwise section also suffers a rotation which to 
a first approximation is the same as that in the preceding case with ribs 
running chordwise. 


For swept-back wings, then, it appears that relative to the centre there is a 
rotation of the chordwise section caused by both torsional and flexural deformations 
of the spar irrespective of whether the ribs run in the chordwise or in the normal 
direction. In the latter case there is an additional camber change in respect of the 
flexural deformation of the spar. In a discussion of control reversal effects on a 
swept-back wing it therefore seems that some differentiation should be made 
according to the particular rib arrangement adopted. 


What evidence there is suggests that quantitatively the effect of the camber 
change is relatively small in the case of the normal rib arrangement. Also, it is 
probably true to say that of the two the chordwise arrangement is the more usual. 
In what follows distortion effects are therefore confined in the first place to rotation 
of the chordwise section and the general effect of the camber change is discussed 
later, as a separate item. The analysis presented in the Appendix, upon which the 
calculated results quoted are based, also excludes the camber effect, being originally 
developed with the chordwise arrangement in mind. The camber effect could of 
course be included. 


2. AILERON REVERSAL EFFECTS 


When aileron is applied the aeroplane rolls as a direct result of the differential 
lift load increments produced on the two halves of the wing by the up- and down- 
going ailerons. Rolling performance, by which is meant the general response in 
roll, is a function of the wing-aileron system, the inertia of the aeroplane about 
the rolling axis, and the forward speed. It is convenient to separate inertia from the 
other effects, which together then constitute the rolling power: a study of rolling 
performance then resolves itself into a study of rolling power, regarded as the 
criterion for aileron action. In other words, rolling power is the effect obtained from 
control deflection in this particular case, the effect that is modified by aeroelastic 
distortion of the wing. Before considering the effects of wing distortion it is necessary 
to establish what shall be taken as a measure of rolling power. 
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Starting from the steady flight condition the application of aileron generates 
a rolling moment and the aeroplane commences to roll, the acceleration in the roll 
being governed by the inertia of the aeroplane. As rolling velocity increases an 
opposing rolling moment proportional to the rate of roll is set up, this being the 
damping in roli associated with the standard stability derivative /,. Eventually a 
steady state is reached when the applied rolling moment is exactly balanced by the 
damping moment: the aeroplane then continues to roll with constant velocity, this 
being the maximum velocity achieved in the roll. The whole process can therefore 
be split up into three stages: the initial stage during which aileron is applied, the 
unsteady stage in which aileron angle is constant but rolling velocity is increasing 
and finally, the steady stage where rolling velocity is constant. 


To make an exact analysis of conditions in the initial stage would be a 
complicated matter as the rate of application of aileron, the rate of generation of 
wing distortion, and the rate of generation of rolling velocity would all have to be 
considered simultaneously. A considerable simplification is effected if it is assumed 
that the aileron angle and the accompanying wing distortion are achieved instan- 
taneously before the aeroplane begins to roll appreciably. The initial stage is then 
in effect compressed into a single instant where the applied rolling moment is the 
result of aileron deflection and wing distortion and is balanced by an inertia moment 
appropriate to the acceleration in roll. Any additional accelerations due to the 
variation of wing distortion with time are assumed negligible in comparison: the 
condition is then effectively a steady state as far as distortion is concerned. Rolling 
velocity being zero, there is no damping moment. Wing distortion is in turn the 
result of the applied aerodynamic loads (due to both aileron and wing distortion) 
and the inertia loads. 


Even with the assumption of instantaneous application of aileron this initial 
condition is probably fairly representative of actual conditions at the beginning of 
the roll. It was, in fact, the practice until fairly recently to simplify the condition 
still further by postulating a static case in which the aeroplane is imagined fixed 
at the centre to prevent rolling, thus deleting the inertia loads. Conditions are 
otherwise the same as in the initial condition described above. 


Recently it has been realised that the static case is perhaps a little too 
unrepresentative of the actual conditions of rolling and emphasis has therefore been 
laid more upon the final steady rolling condition, which has the advantage that it 
can be analysed as an exact dynamic condition without any simplifying assumptions 
and is at the same time much less complicated than the initial or unsteady stage. 
In the steady rolling condition there are no inertia effects, the applied rolling moment 
due to aileron deflection and wing distortion being balanced by the damping moment. 
Wing distortion is in this case the result of the aerodynamic loads due to aileron, 
wing distortion, and damping in roll. 
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One may therefore take as a measure of rolling power either the rolling moment 
in the initial condition (hypothetical but fairly representative), the rolling moment 
in the static condition (hypothetical and less representative), or finally, either the 
rolling moment or rolling velocity in the steady rolling condition (completely 
representative). These three conditions are for convenience summarised in Table I. 


For the avoidance of undue loss of rolling power due to distortion effects it 
has been the practice to ensure that the speed at which reversal occurs, when 
rolling power becomes zero, is well above the maximum design speed of the aero- 
plane. By such means it is guaranteed that over the full speed range there is at 
least positive rolling power. At reversal the three conditions postulated above 
become one and the same. If there is no rolling moment in the initial condition 
there is no acceleration and there are consequently no inertia effects; and similarly 
in the steady rolling condition if the rolling velocity is zero there is no damping. 
Both these conditions are then identical with the static condition with zero rolling 
moment, in which both inertia and damping effects are absent. Rolling moment 
due to aileron deflection is exactly balanced by the rolling moment due to wing 
distortion, and the wing distortion itself is the result of the combined aerodynamic 
loads due to aileron deflection and wing distortion. The reversal conditions may, 
on the same lines as the others, be summarised as in Table II on page 11. 


TABLE I 


Steady rolling | 
condition | 


Initial condition | Static condition 


Zero velocity and | ' 
acceleration | Zero acceleration | 


Rolling motion | Zero velocity 


Applied rolling 
moment due to 


| 


Wing distortion 
due to loads 
from 


Criterion for 
rolling power 


balanced by 


Aileron deflec- 
tion, wing 
distortion, and 
inertia 


Rolling moment 


Restraint at | Damping 
aileron and wing | Inertia moment contee moment 
distortion 


Aileron deflec- 
tion and wing 
| distortion 


Rolling moment 


Aileron deflec- 
tion, wing 
distortion and 
damping 


Rolling moment 
or rolling 
velocity | 


| 
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Having established what is to be taken as a measure of rolling power, the 
effect of wing distortion on the rolling power of a swept-back wing may now be 
considered. Sweepback of course not only introduces additional distortion effects 
but also alters the aerodynamic characteristics of wing and aileron, all of which 
influence the rolling power. The purpose of this paper is not to assess the 
overall effect of sweep on rolling power but merely to provide a clear physical picture 
of the effects produced by distortion on a swept-back wing and of the conditions 
promoting control reversal. The purely aerodynamic aspect of sweep does not 
therefore enter specifically into the discussion. 


TABLE II 


| Reversal condition 


Rolling motion Zero velocity and acceleration 


| Applied rolling moment due 


| to aileron and wing distortion Zero 
Wing distortion due to loads | Aileron deflection and wing 
from | distortion 


| 
Starting with the simplest case of all, that of a swept-back wing in the static 
condition, aileron deflection produces aerodynamic loads which for, say, the port 
half of the wing with down aileron, are as shown in Fig. 5. These “aileron” loads, 
which are distributed over the length of the aileron, are 


A—upward lift loads applied at the aerodynamic centre, 


B—nose-down moment loads applied in the chordwise plane. 


Attention is confined to the port wing, effects on the starboard wing with up aileron 
being simply reversed. 


Assuming for the moment the flexural and aerodynamic axes coincident, these 
aileron loads produce a nose-down twisting and an upward flexing of the wing, both 
of which invoive an aerodynamic wash-out. This introduces the distortion loads C 
pictured in Fig..6, downward or negative lift loads distributed over the whole wing 
and increasing outwards from the centre. In the equilibrium condition the wing 
distortion is then the result of both aileron and distortion loads, as follows :— 
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PORT WING 


| 
AILERON 


LOADS DISTRIBUTED | AERODYNAMIC 


OVER AILERON PORTION. | CENTRE 


| 


CENTRE 


CHOROWISE SECTION 
THROUGH AILERON. 


Fig. 5. 
“Aileron” loads. 


Torsional distortion =nose-down twisting due to component of loads B in normal 
direction. 


Upward flexing= flexing due to component of loads B along spar + flexing due to 
loads A — flexing due to loads C. 


The distortion loads C at the same time apply a rolling moment opposed to that 
produced by the aileron loads A, and in the equilibrium condition 


Net rolling moment =rolling moment due to loads A 
— rolling moment due to loads C. 


Since the aileron loads increase with speed, so also does the wing distortion. 
Loads C, which increase with distortion and speed, therefore increase relative to 
loads A until reversal is reached, when the net rolling moment becomes zero. 


It is interesting to note that the relative amounts of torsional and flexural 
distortion change with speed. At low speeds the distortion loads C are small 
compared with the aileron loads A, whereas at reversal they are equivalent in so far 
as they produce a numerically equal rolling moment. To a rough approximation 
one may therefore say that 


Upward flexing at low speeds=flexing due to component of loads B along spar 
+flexing due to loads A. 


Upward flexing at reversal= flexing due to component ‘of loads B along spar. 


Torsional distortion, on the other hand, is directly proportional to the aileron 
loads B over the whole speed range. 
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PORT WING 


| 


CENTRE 


CHOROWISE TWIST OUE TO 
WING DISTORTION. 


LOADS ‘C’ DISTRIBUTED 
OVER WHOLE WING. 


Fig. 6. 
Distortion loads—aileron. 


It is apparent, therefore, that the ratio of flexural to torsional distortion 
decreases with speed, this being particularly so in the case of moderate sweepback 
where the flexing effect of the component of loads B along the spar is small compared 
with that of loads A. In such a case the upward flexing at low speeds is almost 
entirely due to loads A while the flexing at reversal is practically nil. It is, of course, 
not quite true that at reversal there is no flexural distortion due to lift loads, because 
although loads A and C give equal and opposite rolling moments about the centre, 
the fact that they are not distributed in the same way along the wing means that 
they may still provide some flexural moment at intermediate points. This, however, 
is a secondary effect and would not alter the general conclusions already reached. 


If the flexural and aerodynamic axes are not coincident then in the equilibrium 
condition there will be some torsion about the flexural axis, due to the net upward 
lift loads (A—C) acting at the aerodynamic centre. With a flexural axis in front 
of the aerodynamic axis this torsion results in 2 further nose-down twisting of the 
wing, thus increasing the distortion effects already present. A rearward flexural 
axis would of course have the opposite effect. In any case, the effect of the position 
of the flexural axis, since it depends upon the net lift loads (A—C), must pro- 
gressively decrease with speed and will be small at reversal. One would expect, 
therefore, that at speeds below reversal a forward flexural axis should reduce and a 
rearward flexural axis increase the net rolling power, but that neither should greatly 
affect the reversal speed. 


Finally, there are the additional effects introduced when considering cases other 
than the static condition: namely, rolling inertia loads in the initial condition and 
damping loads in the steady rolling condition. As far as flexural distortion is 
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concerned, the effect is much the same in either case; addition of either inertia or 
damping load results in zero net moment about the longitudinal axis at the centre, and 
flexural moment on the wing is consequently reduced, as in the reversal condition of 
the static case. One may expect to find, therefore, that flexural distortion plays a 
smaller part in the more realistic initial and steady rolling conditions than it does 
in the static condition. Since the damping and inertia loads act downwards, torsional 
distortion will be increased with a flexural axis behind the aerodynamic or inertia 
axis and reduced if it is in front. These departures from the static condition will 
be most apparent at low speeds; as reversal is approached the three conditions, 
as already argued, become one and the same. 


If there is any camber change due to flexural distortion it will, if the wing is 
flexing upwards, be ‘a negative change, involving additional negative lift loads at 
the aerodynamic centre together with nose-up pitching moment loads in the chord- 
wise plane. These additional loads are directly opposed to the aileron loads A and B 
and, in fact, the effect of the negative camber change is mostly simply seen as a 
plain reduction in the positive camber induced by the down-turned aileron: except, 
of course, that the camber change due to flexural distortion occurs over the whole 
wing instead of over the aileron portion alone. The net effect can be represented 
as an effective scaling down of the aileron loads, resulting in a corresponding 
reduction of both distortion effects and rolling power. As the camber change is 
proportional to the flexural distortion its effect will decrease with speed and should 
therefore be small at reversal. 


Aileron reversal effects, that is rolling power and the associated reversal speed, 
can be estimated for a swept-back wing by methods similar to those used for 
unswept wings, the analysis being merely extended along the lines already suggested 
to cover the effects of flexural distortion. Broadly speaking, a choice of two 
methods is offered; either the semi-rigid method in which arbitrary modes of 
distortion are assumed and the absolute displacements at some reference section 
obtained by a direct single-stage solution, or the more exact method in which the 
modes of distortion appropriate to the actual loading are deduced from the elastic 
characteristics as the result of an iterative process. 


An alternative method representing a compromise between these two was 
developed by the author, originally for the estimation of elevon reversal effects, but 
later applied to aileron reversal effects. The method is presented in the Appendix 
to this paper, and a general description of its nature will suffice here. Chordwise 
rotation (relative to the centre) due to both torsional and flexural distortion is 
represented as a series function of the spanwise co-ordinate, the series consisting 
of a linear term plus a Fourier series. Using the Schrenk approximation suitably 
modified to cater approximately for the effect of sweepback, the aerodynamic loads 
are then obtainable as functions of the series coefficients; hence the distortion 
itself is obtained (still as a function of the series coefficients) by an “exact” treatment 
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using the applied loads in combination with the elastic characteristics of the spar. 
The chordwise rotation due to this distortion is then equated to the original series 
function to give an equation (called for convenience the “twist” equation) in the 
series coefficients. 


By applying the “twist” equation at as many stations as there are 
coefficients a direct solution is provided for the coefficients and hence the 
distortion and corresponding reversal effects. The method is in fact analogous to 
the standard Glauert method for determining aerodynamic loading by lifting-line 
theory. With an infinite number of coefficients the method is a true elastic one 
and with some justification may be therefore termed an “elastic wing series solution”; 
but in practice the value of the method obviously depends upon a good approxi- 
mation being obtained with relatively few coefficients. 


This principle is demonstrated by the results of calculations made upon a 
hypothetical V wing. Calculated reversal speeds were obtained using one to five 
terms in the series function and the results are given in Table III. 


The calculations were performed only to slide-rule accuracy, but it is evident 
that in this particular case the two-term solution is sufficiently accurate, while even 
the one-term solution gives quite a fair result, underestimating the reversal speed 
by only 4 per cent. The one-term solution, which involves only the linear term 
in the series function, corresponds to an assumption of linear distribution of the 
chordwise rotation due to distortion: the extent to which the one-term solution is 
in error therefore represents the degree to which this assumption is fallacious. 


In Fig. 7 the calculated values of rolling moment coefficient by the one-term 
solution for the same wing in the static condition are plotted against speed and 
compared with the values given bya semi-rigid method. Aerodynamic characteristics 
were the same for both and the results are therefore directly comparable. In the 
semi-rigid method a linear mode was assumed for the net chordwise rotation due 
to both torsional and flexural distortion; torsion due to lift loads about the 
flexural axis was ignored. 


It is interesting to note that over the first two-thirds of the speed range the 
semi-rigid method provides a good approximation, but the reversal speed is under- 
estimated by 14 per cent. compared with the one-term solution and by 17 per cent. 
compared with the four-term solution. 


TABLE III 
Number of terms 1 2 | 3 4 5 
Reversal speed, m.p.h. 299 314 314 312 | 312 | 
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Fig. 7. 


Static rolling moment for V wing. 


3. ELEVON REVERSAL EFFECTS 


Elevon reversal effects are applicable only to a tailless aeroplane and refer to 
the control reversal effects associated with the use of the elevon as an elevator. 
Because of the fairly recent development that has taken place in tailless aircraft 
the subject is relatively new compared with aileron reversal effects, although the 
essential principles underlying the two phenomena are the same. There are, however, 
certain interesting differences, not the least of which is the relative importance of 
the two effects. Although aileron reversal would be a highly undesirable feature, 
particularly on a military aeroplane, elevon reversal would be completely disastrous 
on any type, involving in its extreme form inability to pull out of a dive. 


When flying at high speed it may not always be absolutely necessary to execute 
lateral manceuvres, but it must invariably be essential to be able to trim the aeroplane 
without danger of structural failure due to large wing distortion consequent upon 
use of the elevon; it is equally important that elevon effectiveness should not 
be reduced to the extent that the available elevon movement is insufficient for trim. 
In short, because of the natural juxtaposition of earth and sky symmetric manceuvres 
in the vertical plane must attain in the limit greater significance than lateral 
manceuvres in the horizontal plane and, apart from manceuvres, the aeroplane must 
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—- 

© DOWN-~ TURNED ELEVON al AERO.CENTRE. 
Z 
CHORODWISE SECTION 

LOADS A + THROUGH ELEVON. 

Ay 
LOADS 8’ 
CHORODWISE SECTION 


INBOARD OF ELEVON. 


Fig. 8. 
“Elevon” loads (no-lift condition). 


be capable of being maintained in the steady flight condition. It is fortunate that 
from the nature of things, as will be seen later, a greater margin is automatically 
guaranteed for reversal of the elevon used as an elevator, than for reversal when 
used as an aileron. 


As rolling power is a measure of aileron effectiveness, so there is an analogous 
elevon power representing the effectiveness of the elevon used as an elevator in 
symmetric manceuvres. Again, as with the aileron, a complete representation of 
elevon power would be based on a dynamic condition, as for instance rate of pitch 
or number of “g” achieved in a steady pull-out. These effects depend ultimately 
upon the power of the elevon as a pitching moment or C,,. producer, and this can 
be represented conveniently by a static condition analogous to that used for the 
rolling power of the aileron. 

Static elevon power is thus defined as the additional C,,. per unit elevon angle 
in the static condition and is represented directly by the total pitching moment 
produced on the wing with deflected elevon when the wing is held at the centre with 
the mean no-lift line of the wing-elevon combination parallel to the air stream. Any 
extraneous symmetric effects due to initial camber or twist of the undistorted wing 
are excluded, only the symmetric effects arising from the elevon deflection 
being considered. 
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CENTRE 


LOADS C 


OUTER PORTION 


Fig. 9. 
Distortion loads—elevon (no-lift condition). 


For this condition the same qualitative analysis of loads and distortions can 
be made as was done for the aileron. The “elevon” loads due to downward elevon 
deflection are shown in Fig. 8 for the no-lift condition, and consist of 


A—lift loads symmetrically distributed on the two halves of the wing, up over 
the elevon portion, down over the inner portion. 


B—nose-down moment loads confined to the elevon portion, applied in the 
chordwise plane. 


Assuming for the moment, as before, that flexural and aerodynamic axes 
coincide, these elevon loads produce nose-down twisting and upward flexing as on 
the port wing in the aileron case. The resulting symmetric aerodynamic wash-out 
provided by both these effects gives rise to the distortion loads C which for the 
no-lift condition are as pictured in Fig. 9: upward lift loads over the inner portion, 
downward lift loads over the outer portion. In the equilibrium condition the wing 
distortion, as before, is finally composed of 


Torsional distortion =nose-down twisting due to component of loads B in normal 
direction. 
Upward flexing =flexing due to component of loads B along spar + flexing 
due to loads A—flexing due to loads C. 
The C,,. produced by the elevon deflection alone is represented by the pitching 


moment due to the elevon lift loads A and moment loads B of Fig. 8. These are 
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additive and give a nose-down or negative C,,.. Distortion loads C of Fig. 9, 
give a positive Cm, and the net C,,. is therefore reduced numerically from what it 
would be for the rigid wing. In the equilibrium condition 


Net Cno=(Cmo due to moment loads B)+(C,,. due to lift loads A)- (Cn. due 
to lift loads C). 


This reduction in power is progressive with speed, as in the case of the aileron, 
and a reversal condition is therefore possible in which the net elevon power 
becomes zero. In the reversal condition interesting differences arise between 
aileron and elevon. 


For the aileron the lift loads C due to distortion balance in effect the lift loads A 
due to control deflection to produce zero rolling moment. To produce zero pitching 
moment in the elevon reversal condition, the lift loads C due to distortion 
must balance both lift and moment loads (A and B) due to control deflection: 
therefore distortion must be greater, involving a higher speed. From qualitative 
considerations alone, therefore, it can be deduced that the elevon reversal speed on 
a given tailless aeroplane must be greater than the aileron reversal speed. 


The second important difference lies in the qualitative effect of flexural distortion. 
Since at reversal the distortion loads C provide a pitching moment numerically 
equal to that given by the elevon lift and moment loads together, then loads C 
must be appreciably greater than loads A. It is therefore quite possible for the 
net downward lift loads (C—A) to be greater in their effect than the component 
of loads B along the spar, with the result that the wing flexes downwards. In the 
case of the aileron the wing was judged to flex upwards, although to a relatively 
decreasing extent, over the whole speed range up to reversal: with the elevon, 
however, the extension of the speed range consequent upon the higher reversal 
speed results in an increasingly downward flexed wing over the upper end of the 
speed range. As the upward flexed wing reduces elevon power, so the downward 
flexed wing augments it and the effect of flexural distortion, therefore, is to delay 
reversal, opposing in fact the effect of torsional distortion. At the low speed end, 
where the wing is flexing upwards, flexural distortion acts in the same sense as 
torsional distortion to reduce elevon power. Since flexural distortion delays reversal, 
it is conceivable that if the flexural stiffness is low enough, compared with the 
torsional stiffness, reversal might be delayed indefinitely. 


The additional effects of flexural axis position and of any camber change due 
to flexing will be similar to what they were in the case of aileron reversal, Except, 
that in so far as they depend upon the net lift loads the effect will be 
reversed at the high speed end. A forward flexural axis position, for example, will 
be bad as it was for the aileron at the low speed end, but will change to a good 
effect at, or near, reversal. 
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Any dynamic conditions involving the addition of “g” effects on the wing, that 
is inertia loads and associated wing lift, if the weights, as is usual, are concentrated 
more towards the centre than the lift, will result in upward flexing of the wing and 
a corresponding wash-out. In the pull-out this will assist the upward-turned elevon. 
On the other hand, curvature of the flight path induces a wash-in which 
opposes the elevon. 


In addition to the question of elevon power there is the problem of trim in 
the steady symmetric flight condition. The two are related in the sense that the 
lower the elevon power, the more elevon angle would be required for trim and, 
proceeding to the limit, an infinite elevon angle would, theoretically, be required 
at the reversal speed. It is advisable, therefore, to ensure that at the maximum 
speed of the aeroplane there is sufficient angle available for trim, with a comfortable 
margin, even though the maximum speed may be well below the reversal speed. 


In the steady symmetric flight condition there are, as in the static condition 
envisaged in the consideration of elevon power, loads due to elevon deflection and 
wing distortion; but in addition there are also loads due to basic (or initial) camber 
and twist of the wing and loads due to gravity and the associated wing lift. Wing 
distortion is the net result of all these loads. For a constant attitude of flight, say 
level flight, loads due to gravity and the associated wing lift remain constant with 
speed and have only a small effect at high speed. Loads due to initial wing camber 
and twist, however, increase with speed and in consequence may have a major 
effect. A wing, for instance, may be given considerable initial positive camber for 
design reasons and in order to cancel the negative pitching moment due to this 
camber an initial wash-out may also be applied, so that at high speed the elevon 
angle required for trim should be small. 


The initial wash-out, in fact, corresponds to tail setting on a tailed aircraft. 
At high speed the lift loads due to this initial wash-out can easily exceed the effect 
of gravity and wing lift and produce a downward flexing of the wing, as opposed to 
the upward flexing which would normally result if there were no initial camber and 
twist. Wash-in induced by the downward flexing of the wing requires a corresponding 
negative or upward increment of elevon angle for balance. The sign of the increment 
of elevon angle required to trim due to wing flexibility, therefore, may be acutely 
dependent upon the degree of initial wing camber and twist present. Particular 
attention is drawn to this feature because in the author’s experience the possibility 
of downward flexing of a tailless wing at high speed seems not always to be realised. 


The elastic wing series solution described in the preceding section is equally 
applicable to the estimation of elevon reversal effects, for which, as already stated, 
it was Originally developed. In the presentation of the method given in the Appendix 
to this paper the two applications are given jointly. The only differences that arise 
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are concerned with the symmetric and asymmetric nature of the two phenomena. 
In the case of elevon reversal effects, for instance, allowance must be made for the 
change in mean wing no-lift line due to the symmetric distribution of wing distortion 
and elevon angle. In the consideration of elevon power attention is confined to the 
static condition. The method is also easily extended to include the effect of loads 
due to initial wing camber and twist, gravity, and wing lift in the steady symmetric 
flight condition. This is also included in the Appendix. In this case the elevon 
angle required for trim enters as a further unknown but is provided for by the 
introduction of an additional equation representing the general pitching moment 
equilibrium of the complete aeroplane. 


Calculations have been made by this method for the same hypothetical wing 
as was used for aileron reversal effects. In the first place convergence tests were 
made, as before, to examine the effect of the number of terms taken in the series 
function. In this case, however, the tests were made by calculating not reversal 
speed, but the increment in elevon angle Ay required to trim due to wing flexibility 
in the level flight condition at 500 m.p.h. E.A.S. The results, given in Table IV, 
indicate that at least three and preferably four terms should be used, although a 
fair result is given by the two-term solution. The error involved in the one-term 
solution is a direct indication of the lack of linearity in the distribution of the 
chordwise rotation due to distortion in this case. 


By systematic variations of certain of the design parameters, the effect of these 
parameters on the elevon reversal characteristics of the hypothetical wing was also 
investigated. The parameters concerned were flexural axis position, aspect ratio, 
angle of sweepback, wing loading, and the torsional and flexural stiffnesses. In 
the case of aspect ratio the wing span was altered in order to keep wing area 
constant. Also, for both aspect ratio and sweepback the elevun percentage chord 
was varied in order to maintain a consistent representative design. For each set 
of variations the elevon reversal speed V, and the increment in elevon angle Ay due 
to wing fiexibility for trim in level flight at 500 m.p.h. were calculated for a range of 
values of the appropriate parameter, using the one-term solution. Although the con- 
vergence tests showed the unreliability of the one-term solution for absolute values, 


TABLE IV 
Number of terms | | | 2 | 3 4 | 5 
An radians - .0406 | -.0213 - 0223 | -.0232 | — .0236 
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it is still reasonable to suppose that it should predict tendencies fairly well. Figs. 10 
to 16 give the results of these calculations as plots of Vz and An against the value 
of the relevant parameter, whose basic value appropriate to the hypothetical wing 
is indicated in each case. 


It will be noted, incidentally, that the elevon reversal speed of the basic wing is 
600 m.p.h. by the one-term solution. Although this may be considerably in error 
it is evident from the convergence test results that since the higher term solutions 
show the wing capable of being trimmed at 500 m.p.h., the true reversal speed must 
therefore be at least greater than 500 m.p.h. This is much higher than the aileron 
reversal speed of 312 m.p.h. for the same wing, thus confirming the qualitative 
conclusions already reached. 


The effects of the design parameters as given in Figs. 10 to 16 are summarised 
in Table V by designating them broadly as either “good” or “bad.” In most cases 
the changes in Vz and An are consistent, one increasing while the other’ decreases, 
but where they conflict precedence is given to the change in A» on the principle 
that more practical significance is attached to the trim condition than to the 
actual reversal condition. 


BASIC 
L 
--2 fe) -2 “4 
FLEXURAL AXIS POSITION BASIC ASPECT RATIO 
BEHINO 4 CHORD. 
Fig. 10. Fig. 11. 
Effect of flexural axis position. Effect of aspect ratio. 
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Effect of sweepback. Effect of wing loading. 
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It is interesting to note that the effect of flexural axis position is as predicted 
qualitatively. For a given flight condition V, is independent of, and Ay is directly 
proportional to, wing loading. Three diagrams are devoted to the effect of varying 
stiffness. In the first, Fig. 14, both torsional and flexural stiffnesses are increased, 
their ratio being kept constant. Reversal speed is proportional to the square root 
of either stiffness and therefore becomes infinite for infinite stiffness and An 
correspondingly zero. In Fig. 15 the flexural stiffness is increased while the torsional 
stiffness is maintained at the basic value. As the reversal speed falls to 500 m.p.h. 
so the value of Ay becomes theoretically infinite. In Fig. 16 the torsional stiffness 
is increased while the flexural stiffness is maintained at the basic value. Here it is 
demonstrated that, as was anticipated, there is a critical value of the ratio of flexural 
to torsional stiffness below which reversal is completely inhibited. 
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Fig. 16. 
Effect of torsional stiffness. 
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TABLE V 


Effect of increasing the parameter 


Parameter 
Ve An General 
| Rearward position of 
| flexural axis Decreases Increases Bad 
| Aspect ratio Decreases Increases Bad 
Sweepback Increases Increases Bad 
| Wing loading Unaffected Increases Bad 
Values of both stiff- : 
| nesses Increases Decreases Good 
| Ratio of torsional to 
| flexural stiffness Increases Decreases Good 
APPENDIX 


ELASTIC WING SERIES SOLUTION FOR CONTROL REVERSAL EFFECTS 


1. Notation 


An effective local dC,,/dz for twist-without-lift loading (see equation (8) ). 
moment of inertia of aeroplane about the longitudinal axis. 

c wing section chord. 

c wing mean chord. 

"0 mean wing lift coefficient. 

Crs section C,, due to lift-without-twist (proportional to C,) 


9 
by Schrenk, C,./C,=0.5 + (7)(S)a - 


Cis twist-without-lift section C;, due to unit elevon or aileron angle. 
for symmetric elevon application: 
Cic= Gm (6. — 6) over elevon portion and = —4,,¢, elsewhere. 


for asymmetric aileron application: 
Over aileron portion 


and =0 elsewhere 
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additional section C, due to rolling, per unit ps/V. 

pitching moment coefficient at zero lift. 

mean pitching moment coefficient at zero lift. 

additional section C,,, due to initial wing camber. 

additional section C,,. due to unit elevon or aileron angle. 

Co Of plain undistorted wing, due to initial camber and initial twist. 
equivalent air speed. 

flexural rigidity of normal section. 

torsional rigidity of normal section. 

acceleration due to gravity. 

total applied rolling moment (whole wing) per unit aileron angle. 
rolling moment derivative due to rolling. 

rolling moment coefficient = L, / qSs. 

value of /,, for the static rigid wing, that is due to aileron deflection only. 


net elevon power in the static condition, that is additional wing C.,, 
per unit elevon angle. 


value of m,, for the static rigid wing, that is due to elevon deflection only. 
moment about a chordwise axis through the flexural centre. 

moment about a normal axis. 

rolling velocity, radians per sec. 

dynamic pressure 4pV”. 

wing semi-span. 

total wing area. 


. moment about a lateral axis through the flexural centre. 


moment about the flexural axis. 
aeroplane forward speed (ft./sec.). 
elevon reversal speed. 
concentrated weight. 


OUTBOARD 


Moments are positive as shown, with the usual convention of the right-hand screw in the 
direction of the arrow. 
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distributed weight per unit span. 
mean span loading = W /2s. 
total aeroplane weight. 


SECTION AERODYNAMIC CENTRE 


spanwise distance from the centre, measured laterally. 

distance outboard along the flexural axis. 

stick fixed static margin. 

distance of section aerodynamic centre behind a fixed lateral axis. 
distance of section aerodynamic centre at x behind a lateral axis through 
the flexural centre at x,, inboard of x. For weights, Y similarly defines 
the c.g. position. 

spanwise co-ordinate x/s. 

angle of incidence, in definition of a... 

local angle of sweepback of the flexural axis. 

increment in elevon angle required to trim due to wing flexibility. 
increase in section no-lift angle due to unit elevon or aileron angle. 
increase in wing no-lift angle due to unit elevon angle. 

elevon or aileron angle relative to wing in chordwise plane, positive 
downward, assumed the same at all sections. 

torsional twist of the normal section, positive nose-up. 

chordwise rotation due to wing distortion relative to section at centre, 
positive nose-up. 

increase in wing no-lift angle due to wing distortion. 

initial chordwise rotation of the undistorted wing relative to section at 
centre, positive nose-up. 

wing no-lift angle relative to section at centre, positive nose-up, for 
the undistorted wing. 

air density. 

change of slope due to flexure, measured along the flexural axis. 
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2. Chordwise rotation due to wing distortion 
Equilibrium conditions at any section give 
M’=M cos y - T sin y 
T’=T cos y+M sin y. 


The net increment in chordwise rotation due to torsional and flexural 
distortion is 


dé=d cos y — dw sin y. 


Also, dé’= (=) ax 
and ax. 
Substituting for d6’ and dy, and for M’ and T’, and writing dx’=dx/cos y 
gives finally 


dé 
a =K.T+K,M , . . . (1) 


where K,= (4 + EI COS 7 


1 


Equation (1) is the general equation for chordwise rotation due to wing 
distortion and is applicable to all reversal effects. The effective flexibilities K, 
and K,, will vary in general from section to section. 


3. Elevon reversal effects on a tatlless aeroplane 


The effects considered are symmetric and account must be taken of the position 
of the mean wing no-lift line. 


Chordwise rotation due to wing distortion is represented by the series function 
(n=0, 1, 2, 3, etc.) 


where, in the case of n=0, sin nzz is replaced by z. 
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The net incidence of the wing relative to the mean no-lift line is (4, - 8.) due 
to initial twist and (6-4) due to distortion. By the Schrenk approximation 


1 

1 

0 


1 
and ,¢, | sin (c/c) dz. 
0 


The contribution to wing C,,. provided by the distortion is 
| an (6 (c/c) (y / c) dz 
0 
= > a, aks 


1 
where ,¢.= | Am (ng, — sin nzz) (c/ (y/c) dz. 
0 


Adding C,.. for the plain wing (including the effect of the initial twist @,) and 
m,xn for elevon deflection then gives the total wing C,,... The condition for pitching 
moment equilibrium is therefore 


referred to generally as the “trim” equation. 


Wing distortion effects are represented by applying equation (1) with the 
appropriate conditions. The net component torque and bending moment due to 
aerodynamic loads and weights are, at station x,, 


T 3 Cmoqe*dx - CyrqcYdx + waYdx+ (5) 
M= (x x,)dx- Wa (x-X,) dx- = (x-x,). (6) 


1 1 


The total section C,,. and section C;, due to all the various effects, including wing 
distortion, are 


Cy = Crs + an (6, 6,)+ an (6 6) (8) 
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Substituting (5), (6), (7) and (8) in equation (1) gives finally 


dé = 


where, for station x,=Z,5, 


A,=K. { Cuow ~ 2) dz + 


1 


Ky, Am (8 9.) =") dz 


B,=K, { - (%) a} 


1 


fi 
nhs = an (6, — sin nzZ) (c/c) (Y /c) dz- 


- Ky | Am — sin (c/ Cc) (2) az 


1 


K, and K,, are the effective flexibilities at station x,. 
Integrating equation (9) gives the chordwise rotation due to wing distortion as 


6=qs*C? (A, + Zan. 


where A,= A,dz, etc. 


0 


and ,(,= 
0 
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This must agree with equation (2), and equating the two gives finally the “twist” 
equation as 


Dan sin nzz } +A,+B,C,+C.=0. ; 


For a given speed, g and C, are known. The “twist” equation (11) is applied at 
as many stations as there are coefficients a,. With the “trim” equation (4) thic 
provides a set of simultaneous linear equations for the solution of the coefficients a, 
and elevon angle 7. 


Net elevon power in static condition 
The net elevon power (C,,. per unit elevon angle) due to elevon deflection and 
wing distortion is 


where the coefficients a, give the chordwise rotation due to distortion arising from 
application of unit elevon angle in the static condition. 


The “twist” equation is equation (11) with the A and B functions omitted to 
exclude all effects other than those due to elevon deflection and wing distortion. 
With »=1, the “twist” equation then becomes 


Sa, { (=) sin naz } +C,=0. 


For a given speed the “twist” equation (13) provides a set of simultaneous 
linear equations for the solution of the coefficients a,, which when substituted in 
equation (12) give the net elevon power m,,. 


The elevon reversal condition is given by putting m, in equation (12) equal 
to zero. With the “twist” equation (13) this provides a solution for the coefficients 
a, and the speed V = /(2q/p) (see later under Aileron reversal condition). 


4. Atleron reversal effects 

These effects are asymmetric and there is therefore no change in the mean wing 
no-lift line. Consider, for example, the port wing with unit down aileron angle. 
Static condition 


The only loads which arise are those due to aileron deflection and wing 
distortion. 
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The total rolling moment for the whole wing due to aileron deflection alone is 


2 Cy -qcexdx =1,2qSs, 
0 


where len = | Cie (¢/c) zdz. 
0 
The total rolling moment due to chordwise rotation is similarly 


2 | GmOqcxdx = Qn 
0 


1 
where now | sin nxz (c/c) zdz. 
0 
The net applied rolling moment is then L, =1,qSs, 


The “twist” equation is exactly the same as equation (13) for elevon power, except 
that the functions are for asymmetric aileron effects instead of symmetric elevon 
effects. 9 being now zero, ,¢, is derived with ,¢,=0 and C, is also derived using 
the appropriate values for C;.. 


For a given speed the “twist” equation provides a set of simultaneous linear 
equations for the solution of the coefficients a,, which when substituted in equation 
(14) gives the net rolling moment coefficient /,. 


Initial condition 


Additional downward loads arise due to inertia appropriate to the rolling 
acceleration dp/dt, the inertia load on a weight w being (w/g)xdp/dt. Also, since 
dp/dt=L,g/A=qSsl,g/A, the inertia loads can be expressed in terms of the net 
rolling moment coefficient /,,. 


By a process similar to that followed for the other effects, the effect of the 
inertia loads on the wing distortion is given by an additional term in the “twist” 
equation, which becomes 


where E,= E dz, 
0 


A 


32 


). 


CONTROL REVERSAL EFFECTS ON SWEPT-BACK WINGS 


For a given speed the rolling moment equation (14) and “twist” equation (15) 
provide a set of simultaneous linear equations for the solution of the coefficients a, 
and the net rolling moment coefficient /,. 


Steady rolling condition 


In addition to the loads due to aileron deflection and wing distortion there are 
also aerodynamic loads due to damping in roll. These loads are proportional to 
the rolling velocity p and can be represented by an additional section C, given by 


AY 
The “twist” equation now becomes 


{ ake - sin } +C,+D, ==(), . (16) 
where D, = j, D, dz, 
0 


ri 
and D,=-K.| Cin Ku Crp (c/0) az. 


1 


Rolling moment due to damping is given by the standard stability derivative /, 
as equal to gS.2s./,.ps/V. For equilibrium of rolling moments the rolling 
moment equation then becomes 


For a given speed the rolling moment equation (17) and “twist” equation (16) 
provide a set of simultaneous linear equations for the solution of the coefficients a, 
and the rolling velocity p. 


Aileron reversal condition 


At reversal the static, initial and steady rolling conditions become identical. 
Putting /,=0 in equation (14) or p=0 in equation (17) gives the same rolling 
moment equation, 


Similarly, putting /,=0 in equation (15) or p=O in equation (16) gives the same 
“twist” equation as for the static condition, namely 
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4\. 
Xa, (5) SIN } +C,=0. (19) 


The rolling moment and “twist” equations provide a solution for the coefficients 
a, and the speed V = /(2q/p). As in the case of elevon reversal the equations are 
non-linear in a, and gq, but the coefficients a, can be eliminated to give a direct 
solution for g. For aileron reversal, for instance, the solution is given by 


oS2 1$2 mS2 lie 
a a a 
(20) 
m+1%o m+1%1 m+1%m 


where, in general, 


4 
sin s7z at station 7=Z,, 
and ,C, is the value of C, at station z=z,. (m, rand s denote here any integer). 


If this direct form of solution proves too cumbersome it may be simpler to 
calculate the net power for the static condition over a range of speeds and hence 
find by graphical interpolation at what speed the power is zero. 
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CALCULATION OF DOWNWASH 
BEHIND A SUPERSONIC WING 


by 
G. N. WARD 


(Department of Mathematics, The University, Manchester) 


Once the linearised lift problem for a supersonic wing has been solved and the 
disturbance velocity potential or the pressure is known at every point on the 
surface, the complete flow round the wing and the downwash behind the wing may 
be found by the following procedure. 


Let x, y, z be a right-handed system of Cartesian co-ordinates such that x is 
measured in the direction of the undisturbed stream (velocity U) upstream from the 
wing. Then the disturbance velocity potential ¢, such that the disturbance velocity 
u is given by 


u=V¢9, ; ; (1) 
satisfies the well-known equation 


d%@ 


where B?=M?-1 
and M is the Mach number of the undisturbed stream. 


Consider now a wing of infinitesimal thickness at a small finite incidence to the 
stream and lying very nearly in the plane z=0, then it follows from the antisymmetry 
of the system that the potential and its x and y derivatives will be odd functions of z, 
and d9/0z will be an even function of z. Since d¢/0x (the pressure) must be con- 
tinuous except through the wing, it follows that (0¢/%x),_,=0 except over the wing. 
There will be a vortex sheet extending downstream from the trailing edge of the 
wing and lying approximately in the plane z=9, for which the vorticity must be 
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constant along streamlines, which in this case are approximately straight lines 
parallel to the x-axis. Hence the difference of potential between the two sides of 
the vortex sheet must be a function of y only. Also the potential must be continuous 
everywhere except through the wing and the vortex sheet, so that the potential on, 
say, the upper side of the vortex sheet must be equal to the potential at the trailing 
edge of the wing on its upper side (which is assumed to be known) on any line 
y=constant. In the plane z=0, at points not on the wing or the vortex sheet, we 
must have ¢=0 since the potential is continuous at these points. Thus the potential 
is known everywhere on both sides of the plane z=0. 


From the anti-symmetry of the problem we can confine our attention to the 
half-space z >0. The potential ¢ (x, y, z) satisfying equation (2) and having the 
value ¢ (x, y, +0) on z= +0 is given by the known expression due to Hadamard,"’? 


th 1 ) 
(x, y, 2) (x’, y’, +0) x) -B dx'dy’ (3) 


the integration being taken over the area of the plane z=0 for which 
x’ <x-B{(y-y’h +27}! . ‘ (4) 


and where * denotes that the “finite part” of the integral must be taken. 


In the case under consideration ¢ (x, y, +0)=0 except on the wing and the 
vortex sheet, so that the integral is taken over those parts of the areas of the wing 
and the vortex sheet for which equation (4) holds. 


By integrating equation (3) by parts with respect to x’ and remembering that the 
potential is continuous and is zero on the leading edge of the wing, we have 
1 y’, +0) 
(5) 


the integration now being taken over that part of the wing alone for which equation 
(4) holds, since 6¢/0x=0 on the vortex sheet. This gives everywhere in z > 0 and 
the velocities may be found by differentiation. The quantity usually required will 
be (9¢/9z)._,, giving the downwash in the plane of the wing. 


In order to illustrate the method it will be applied to calculate the downwash 
far behind a symmetrical flat delta wing at incidence, whose leading edges lie inside 
the Mach cone from its vertex. The result for this case has been given previously 
by Robinson and Hunter-Tod: *) following some rather complicated analysis. To 
facilitate comparison of the results their notation is used. 


36 


| 
|_| 


DOWNWASH BEHIND A SUPERSONIC WING 


Let y=semi-vertex angle of the wing, »=tan y, A= Bn, 
c=maximum chord of the wing, 
d=distance behind the trailing edge of the wing, 
a= incidence of the wing. 

Take the origin of co-ordinates at the vertex of the wing. 


For this delta wing it is known) that 


(* 2=0 E’ (A) (n?x? — y?)! (6) 


When x > c(1+A)+B\z{|,ie.d>Ac+B|z| , equations (5) and (6) give 


(x, Z)= EF’ (a) (x x’)x'dz’ x 


nx dy’ 
= z(x- x’)x’dx’ x 


#/2 dé 


(y—nx’ sin 6) +z?) [ B? (y - nx’ sin 6 — 


x 


The evaluation of (0¢/z),_, from this integral is a little tedious in the general 
case, but for y=z=0 the integral assumes a simple form; in this case, the integral 
with respect to @ becomes 


o sin?@ [ (x x’? — B’z? — sin 76) 


As z-0, this integral is of the form 


and in order to calculate (0¢/9z)._, we require only the term A,, because of the 
factor z in (7). Now 


_ lim _22_ (*” do 
A_, n*x/2 sin? 6+ 2? . . . (10) 


Therefore, if K=Ax’/(x - x’), 


, ( l ) 
\, sin? 6(1-k? sin? 6)! sin? 6 


(" 2 sin* 


= K (k)- E (k). 
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Hence 
(x, 0,0) _ 22Un* K (k)- E(k) 


0z rE’ (a) dx 


Ac 
_ ATK (k)-E(K)) 
A 


/ 


which is the result given in Ref. 3 for d>Ac. The corresponding value for 
0<d< Ac is easily obtained by this method and is given in Ref. 3. 
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ESTIMATION OF THE EFFECT OF 
A PARAMETER CHANGE ON THE 
ROOTS OF: STABILITY EQUATIONS 


by 


K. MITCHELL, Ph.D., A.F.R.Ae.S. 


(Lecturer in Mathematics at King’s College, Newcastle upon Tyne) 


SUMMARY 


A method is given for calculating approximately the changes in the roots of a 
stability secular equation caused by a change in any of the parameters involved. 
General formule are given applicable to any quartic equation, and special formule 
are also given applicable to the stability of an aeroplane: lateral stability in the 
text, and longitudinal stability in an appendix. 


The method of using the formule is illustrated by applying them to a particular 
calculation of the lateral stability of an aeroplane, and a check of the results is made 
by comparing the predicted approximate changes with those calculated by solution 
of the modifiéd period equations. It is shown that the formule are reliable, for this 
typical case, for any reasonable changes in any parameter other than n,. 


If the changes in the derivatives are made equal to the probable error with 
which they can be measured, the formule enable us to evaluate the probable errors 
of the roots. These are found to be considerable, and to arise mainly from uncer- 
tainties in y,, n,, and n,: if these could be reduced to 0.03 in y, and 0.006 in the 
others, the uncertainties in the roots would be reduced to some ten per cent. of their 
values, except for a larger uncertainty in the root corresponding to the 
slow spiral motion. 


1. INTRODUCTION 


There can be few people who have completed any calculation in the field 
of aeroplane lateral stability without feeling some doubt as to the physical meaning 
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to be attached to the results. The difficulty resides in our imperfect knowledge of 
the values of the derivatives and other parameters entering the lateral stability 
equations. Some of these are commonly measured in wind tunnels, and may be 
taken as known within the accuracy permitted by experimental technique, and by 
scale effect from tunnel to full scale. For others we use theoretical estimates, based 
on a theory sketchily confirmed by a handful of experimental results. Of other 
derivatives we know littlé, and assume less: we take them to be zero! What 
credence, then, are we to give to the results of our calculations? 


Our doubts could be resolved by better knowledge of the lateral stability 
derivatives, but this is likely to come only slowly, in view of the difficulty of the 
subject, both theoretically and practically. There is, however, another way. We 
can carry out a numerical examination of lateral stability, to see how great in fact 
are the uncertainties in the roots corresponding to the uncertainties of the derivatives. 
Such a survey, carried out by exact solution of the secular equations for lateral 
stability, would be an overwhelming task. Approximate methods must be used, 
but not approximate formule for the roots, such as those given by Bryant and 
Hopwood"), Gandy”), and Mitchell, Thorpe, and Frayn®). These formule 
represent compromises which attempt to give the values of the roots over a wide 
range of conditions, and sacrifice accuracy to range of applicability. What we need 
for the proposed survey is a method of estimating the rates of change of the roots 
with respect to the parameters concerned. One possible method* is found by 
generalising an idea due to Cowley and Skan“). 


Cowley and Skan examined the rates of change of the reots of a polynomial 
equation with respect to its constant term, and went on to develop a method of 
solution of polynomial equations depending upon the integration of the formule 
obtained. The author has never been able to regard this transformation of a 
polynomial equation into a differential equation otherwise than as a retrograde step, 
as a method of solution. We may use it, however, to estimate the approximate 
changes in stability due to a small, but not infinitesimal, change in the constant 
term. In this form, and generalised to include effects due to variations of any 
coefficient, and by a further extension to include those due to variations of any 
parameter upon which the coefficients depend, it is exactly what is desired. 


This method of calculating the rates of change of the stability roots with respect 
to the stability derivatives enables us to solve our original problem, and some others 
of rather greater practical importance. The designer of aeroplanes is often faced 
with the problem of estimating the effect upon stability of some proposed design 


* A similar method is applicable to many types ci stability problems, and the flutter problem 
has already been studied in this way by Buxton and Minhinnick(*); the author was 
unaware of this paper when the above introduction was written. 
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change, or with the converse problem of amending a proposed design so as to bring 
about a desired change in stability. The methods to be developed in this paper 
enable him to solve both problems in a rapid approximate manner and, if it is 
desired to go further and evaluate the changes accurately, the modified coefficients 
in the modified secular equation can be obtained exactly from formule 
to be given later. 


It is thought that this problem is of sufficient importance to the designer to 
justify the preparation of a paper mainly devoted to it. Our original problem, the 
investigation of the probable errors of stability roots, will involve much numerical 
computation: to await the completion of this would delay publication, and we 
return to it only briefly 12 Section 6. It will be necessary to evaluate the rates of 
change of the stability roots with the cerivatives over a large range of possible cases, 
and it is to be hoped that, while answering the somewhat theoretical question 
proposed, the survey of lateral stability that is involved will yield results of more 
direct practical value. Some of the particular formule for the rates of cuange may 
be too complicated to be of other than numerical use: others, it is to be hoped, 
may lend themselves to simplification and approximation, so contributing to our 
general knowledge of the subject. 


It may be added that it is a necessary condition for the validity of the formule 
to be developed below, that they are not to be applied in cases where the stability 
equation has nearly equal roots. It is precisely this condition which has led us to 
confine our attention to lateral stability in the text of this paper, though some of the 
corresponding formule for longitudinal stability are given in an Appendix. In 
lateral stability* equal roots appear only for negative values of the parameter n,, 
or for very small n, coupled with positive /,: both are conditions which should be, 
and usually are, avoided in practice. Our formule are therefore of general usefulness 
in the lateral stability field. In longitudinal stability, on the other hand, the region 
of the stability field in which equal roots occur can hardly be avoided in practice, 
though its stability characteristics are far from perfect. In this case an approximate 
theory already existst, which, in the circumstances, is probably of more use than 
the methods of this paper. The appendix on longitudinal stability, therefore, con- 
tains merely the basic formule for the changes in the coefficients due to changes in 
the derivatives. These, although of little theoretical interest, will be of assistance 
in the setting up of the modified stability equation corresponding to any proposed 
change in design. 


* The only published calculation of the “locus of equal roots,” the curve in the /y—nv 
diagram along which the roots are equal, was given by Mitchell,.Thoroe and Frayn(). 
A few other calculations have been made by the same authors (unpublished), and the 
statement in the text is true as far as present knowledge goes. 


+ See Neumark(®), especially Appendix II. 
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2 VARIATIONAL FORMULA FOR A QUARTIC 


2.1. General 


We now consider the rates of change of the roots of a polynomial equation with 
respect to a parameter upon which its coefficients depend. We take the quartic 
equation 


f (p)=Ap*+ Bp’ +Cp?+Dp+E=0, . 
whose coefficients A, B, C, D, E, are functions of a parameter «, and take standard 
values A,., B., C,, D., E., when this parameter takes a standard value «,. If p 
denotes the value of a root of (2.1.1) when « has its standard value «,, and if p,+dp 


is the value of the corresponding root when «- «,+dx, we have, to the first order 
of small quantities, 


{ A, +(0A /dx)dk } (pot +4p.°dp)+ { B,+(OB/dx) dk +3p.2dp) + 
+ { C, + dx) dk } (p.? + 2p.dp) + { dk } (Do + dp) 
+E, +(0E/d«) dx =0. 


Subtracting the terms of zero order, which vanish separately, we obtain a relation 
between dp and dx. In the case of a real root, p,= - q, this takes the form 


0A 4 oB 3 aC 2 oD + dE 
dq= - dk — —— (2.1.2) 
4A.q° + 2C.g—D, 
For a complex pair, p,= —r+is, we obtain 
dr ids = — dk — , (2.1.3) 


4A, (r — is)’ — 3B, (r- is) +2C, (r-is)—D, 
which can be split into its real and imaginary parts to give formule for dr and ds 
separately. 


These formule express the rates of change of the roots as linear functions of 
the quantities (1/A,)A/9k, etc., the coefficients being functions of the roots alone.* 


* For, if —q, be any root of (2.1.1.), and —q,, —q,, —q,, the remaining roots, for the standard 
case, then 
4A oq,° —3Boq,” + 2C09, — Do = — — — 
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The formule are conveniently written in the forms 


dq 1 B ac (2.1.4) 


K K K 


for a real root —q, with coefficients K,,, etc., and for the real and complex parts 
of the complex root —r+is, with coefficients K,,, etc., and K,,, etc. 


The values of the twenty coefficients K depend upon the values of the roots, 
and upon whether they are real or complex. The cases of four real, two rea! and 
two complex, and four complex, roots therefore require separate consideration. 


2.2. All roots real 
The roots may be denoted by p= — qg;; i= 1, 2, 3,4. We then have 


with similar expressions for the other coefficients K,;:. Also, in general, 
Calculated values may be checked by the relations 


4 


4 4 
= = = Kyic=9 
1 1 


2.3. Two roots real and two roots complex 
The roots will be denoted by —g,, —gq., real, and —r+is, the complex pair. 

We then have 


(2.3.1) 
— 9.) { + } ] 


the remaining coefficients for the real roots being obtained from (2.2.2) as above. 
For the complex pair, we have 

(2r- - q2)/[2 +8 +87} (2.3.2) 


43 


a 
p 
lp 
_| 
4 | 
i=1 | 
i=1 
f 
| 


K. MITCHELL 


and the remaining coefficients may be calculated from the relations 


Kp = 


, (2.3.3) 
etc. Checking relations are 
Kas + 2K +x = Kap + Koop =Kare + 2K | 
Koon +2K-3= { (2.3.4) 
J 


Koos +2K..= — B,/ 
These check by implication all K, up to K,g: to check K,, we use the additional 
identity 


2.4. All roots complex 
The roots are -r,+is,, —r.+is,. Then we have 


(r,—r.)/[ { (r,—r.¥ +(s, }{ (r,—r. +(s,+5.% }] | 


Kye= {(r,- }/[2s,{ + (8, 12)? + (8, +52) }]| 


with similar expressions for K,.. and K,... The remaining coefficients K,,p etc. are 
calculated by formule of type (2.3.3), and the checking relations become 
Krop= K,,0+ Kroo=0 } 
- 4B,/ A. 


checking all save K,,. To check these, we add the identity 


— Koy, =4 — B,”)/ . . (2.4.3) 


2.5. Derivation of checking relations 
Using —q;; i=1, 2, 3, 4; to denote the roots of (2.1.1), we have the partial 
fraction expansion 
A 1 l 
f (p) (4: 492) (91-93) (91-93) 
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Expanding both sides in ascending powers of 1/p, and on coefficients, we 
conclude that 


n 


(q, — (q,- qs) (Gi qs) 


= 


is equal to the coefficient of (- 1/p)'** in the expansion of A./f(p). Thus 


0, >= l | 
Ao, o,=(B,?—A.C,)/ A.’ | 


These lead immediately to the previous checking relations. 


3. APPLICATION TO LATERAL STABILITY 
3.1. Notation 


For the lateral stability period equation we have 


le I, Ig To 


C (* Me Me, in ke 


I, iy is io ic ic iy 


is Io 
D=-y, 


I, Io ic ic 


Cc 


(3.1.1) 


ic Io Ty ic Io i, 


The notation used is that of Bryant and Gates‘), with the replacement of 
their 6, by y., and with the addition of the dimensionless product of inertia coefficient 
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and of extra dimensionless derivatives y, and y,, defined in terms of their natural 
counterparts Y, and Y, by the equations 


Here S is wing area (ft.*), b is span (ft.), U. is forward speed (ft./sec.), p is air density 
(slugs /ft.*), m is mass (slugs) and E the xz-product of inertia (slugs. ft.’). 

Care should be taken not to confuse E of (3.1.2) with the coefficient E of (3.1.1). 


We shall also use the supplementary notation developed by the author‘), 
with minor modifications—the replacement of 2 by 3 in /, and n., now /, and n,, 
and of y, by y.. The symbolism found useful here is as follows: 


L= I,= /is 
Ny /ics n,= N,/ic 
(3.1.4) 
€,= — iz/is, — 
and 
ln, +l1,n, 
1 €,€, 
(3.1.5) 
WHET 


These symbols allow us to express the coefficients of the period equation in the 
simple and convenient forms 
B/A=T,+y, 
C/A=T,+T y.+ (l—ye/ m2) +1 
(3.1.6) 
D/A=V,y, +11, + k (IL, tan ye) 


E/A=k (II, - II, tan y.) 
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j 


3.2. Formule for variations of coefficients 


In calculating changes of coefficients, or the resulting changes in roots, it is 
found convenient to introduce an intermediate parameter «’ and to regard the 
quantities (1/A,)0A/0x, etc., as the products of the quantities (1/A,.) dx’ /d« and 
0A/dx’. Formule for these quantities are given in Table I. The parameters taken 
ALC Yoo Vr Lvs Up, Le Mvs Mr Mo, tan ye L/is, 1/ic ig : the coefficients A, B, 
C, D, E, are linear functions of each of these, except iz, and the table can therefore 
be used to calculate exactly the changes in the coefficients due to any changes, small 
or large, in these parameters. If « is changed by Ax, we deduce that the coefficients 
become A, + 4A, etc., where 


AA=A, (4 54) de. 


To render possible the same accurate calculation of changes in coefficients due 
to changes in i,, we have added two extra rows to the table, giving the rates of 
change with respect to i, and (i,)’, regarded as independent. The changes are then 
calculated by means of 


LA= A. (+ dk ) AS , iy Aig + Aol (3.2.2) 


and similar relations. 


4. NUMERICAL EXAMPLE OF THE USE OF TABLE I 
4.1. Data 


We proceed now to give a numerical example showing the use that may be 
made of Table I. We shall suppose that an investigation is being made into the 
lateral stability of a lightly loaded aeroplane of large wing span flying level at a 
speed near its maximum. Values of the derivatives and other parameters are taken 
to be 

W= -02, y,=y-=0 

l,= —0.06, - 0.42, 1,=0.06 
n,=0.024, n,= 0.03, — 0.048 (4.1.1) 
0.10, 
w,=10, C,=0.2, tany.=0, k=0.1 
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From this the roots of the period equation are evaluated; there are two real roots 
and two complex roots, —q,, —q,, —r+is, where 


q:=0.010 946, g,=4.272 817 
(4.1.2) 


r=0.218 118, s=1.415 952 
In the process of setting up and solving the lateral stability period equation certain 
other quantities are evaluated: these are 
L=42, i,=-06, a,=02, 2n,=0.32 
y,=02, WV=16, ¢,=6,=0 
=4.52, T,=1.464 
(4.1.3) 
B,=4.72, C,=3.968, D,=8.8128 
E,=0.096 


4.2. Calculation of variations of roots 


The suggested application of the formule of this paper is to the estimation of 
the stability change due to a design change, which may affect only a single parameter. 
To make any such calculation, it is necessary to evaluate 4 (n+ 1) of the 20 constants 
K of Section 2.1, where n is the degree of the term of highest degree in the period 
equation in which the parameter to be varied occurs. Here, for completeness, we 
evaluate the complete set at once. Using (2.3.1) and (2.2.2) for the real roots, and 
(2.3.2) and (2.3.3) for the complex roots, we find the values given in Table II. 


TABLE II 


Values of the quantities K, etc. 


K K, K Ky Ky 


Ko 0 | O 0.000 014 | —0.001 254 0.114 579 
K,, » 4.240 006 0.992 321 | -0.232 240 0.054 353 -0.012 721 
K - 0.239 999 0.003 841 | 0.116 114} -0.026 550 —0.050 929 
K, 0.031 403 -0.168 905 | 0.020 599 0.077 914 | -0.026 596 
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TABLE I 
Values of quantities (1 / A.) dx’ /d«,0A/d«’, etc., needed in ¢ chan 


(a) Parameters entering linearly. Coefficient change (exact) AA=A Ax, etc. 


1 dk’ 
YA / dK B/d« NC / dk 
Ye 0 1 
Yo t1/p, 0 0 II 
Yr —1/p, 0 0 II, 
Me 0 0 Yo (i =) (1 
i, — €,€, U 
| 
l 0 n,+ 
iy (1 — sty 
at y: Yp 
) 
n pase 0 1 L+y 
lo € ) 
Ms +1/p, 0 0 II, 
| tan y. —k 0 0 0 
} P2 
€,€ 
(b) Parameters entering other than linearly. Coefficient change (exact), 4A =A, ( + =H 
2 | 1 | 
€,€) 0 Ye (I. + n,)— 1, Yr) +MyYp 
l 
| 1 Ya 0 


(1 — 


(c) As (6). Approximate root change, or approximate’ coefficient change, as (a) above. 


| 1 | | 
inde | | Yo (le (us —Ye)+ Mp | 


| 


ating changes in coefficients and stability roots. 


dk’ 
A, dk 


¢, etc. Root change (approximate) 4g= 


OK’ 


aD /d«’ 


Yr n,+ +k (1—e¢ tan ye) 


> y 
(1- ) 


y 
A 


(1- 721, +k 


hyt (1-%) 


Yo « 


Iho 


(il, — I, tan y) 


Il, 
Do +k NV tan 
D.-k& 


r\A, 0K 


1 dA =) (4 
| 


bok (n, +1, tan ye) 


pok (n, +1, tan ye) 


dE/dx’ 
0 
0 
0 


k (n,—n, tan ye) 
—kW tan ye 
-k (l, +1, tan y.) 
tan ye 


k (II, — IT, tan ye) 


0 


§ 
r, | 
Il, 
I, 
| 
TI, | 
E. 
| 
| 
| | | 
0 0 | 
above. 
| 
| | 0 
Yo | | 
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These satisfy the checking relations (2.3.4) and (2.3.5), the residuals, starting with 
K., and ending with (2.3.5), being 0, - 1, +2, +3, —4, +5, in units of the sixth 
decimal place. 


The next step in examining the effect of a single design change would be to 
evaluate the appropriate quantities (1 / A.) A /0x«’, etc., from a single row of Table I, 
and by multiplying these by the appropriate sets of constants K from Table II, to 
form the rates of change of the roots with respect to the varied parameter. We 
shall multiply the labour involved in so doing 14-fold, by carrying out this calculation 
for variations of each of the parameters. We therefore proceed to construct a 
numerical version of Table I, which we give as Table III. Making only a single 
parameter change, of course, we should need only a single row of this table. 


TABLE III 


Typical numerical version of Table I 


1a | | aD 


| | 
0 1 | 452 | 1464 | 0 
| oy, | +01 
| | | o | 1 | 792 | o | 
, | -t0 | o | 1 | 052 | 1664 | 0 
l, 0.2 0.04 | -0.16 | 
| | +200/3| 0 | 0 42 | -0.06 | 
m | -20/3 | 0 1 | 44 | 084 | 06 | 
wn, | +01 | 0 0 1.6 852 | 0.096 
| tany | ~0.1 0 0 0 16 | 792 | 
Hix | +001 0 | 42 | 2304 | 88128} 0,096 
lic | +0.15 0 0.32| 3.128 | 8.2128 0.096 
i, | -200/3|} 0 0.03 0.606 | 0.024 | 0 


From this we construct the gradients of the roots with respect to the parameters, 
obtained by forming the sums by rows of products of quantities in corresponding 
columns of Tables II and III, and multiplying the final result by the appropriate 
quantity (1/A,)9%«’/d«. For example, taking the top rows of both tables, we find 
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1) { Ox 0+1x0+4.52 x 0.000 014+ 1.464 x (—0.001 254)+ 


+0x 0.114 579 } =0.001 773. 


The quantities so obtained are given in Table IV, which contains also the 
residuals, in units of the sixth decimal place, calculated from the checking formule 


or 1 OB B, 0A 


dg. _ ds B, A.C.—B,? 2A 


TABLE IV 
Table of rates of change of roots with respect to parameters 
Residuals | 


| | | | | 
| | gr | OK | or / OK | / OK | (4.2.1) | (4.2.2) 


yy | 0.001 773 | - 0.022 169 | —0.489 807 | —0.038 269; -—10 | +8 | 
yo | —0.000 112 | 0.134 126 | 0.067 120} 0.019 839 +2 | 


| | | 

| ye | 0,000 991 | --0.005 889} 0.002 -0.065 0 | 


| 


| -0.329 033 | 1.589 883} 0.959 473 | —2.252 313| +30 | 


| 


l, | 0.020 794 | -9.619 996 | -0.200 411} 0.285 446| -24 | +163 | 


1, | 0.183 800 | -0.422 385) 0.303 094) 0.114917| +3] -3 | 
n, | -0.808 503 | —0.212 943 | 0.510 649 | 23.295 571 | —148 | —107 
n, | 0.051 107 | - 1.288 643 0.618 784 | —3.261 287| +32 | +11 | 


n, |—0.451 704 | 0.056 593 | - 3.079 221| 0.191 861| -72 | -21 
| | 0.000 034} 0.009 028 | —0.004 531} 0.069 423) 0 
‘tan y. |-0.090 546 | 0.001 379| 0.044 584| 0.008 598 | +1 
| |-0.000 002 | 0.411 045 | 0.004 479 | 0.002 215) +1) +14 | 
| | 0.000 112} 0.005 440| 0.021 225} 0.097258; +2) 42 

| ip | 0.001 441 | 7.310 889 | —4.656 207 | -0.619 052| -84 | | 
| | | | | 


5. TEST OF THE FORMULA BY CALCULATION OF EXACT AND 
APPROXIMATE CHANGES IN ROOTS 
The entries in Table IV are the rates of change of the roots with respect to the 
various parameters, and to check the validity of the method, we must calculate the 
exact changes in the roots of the secular equation caused by various hypothetical 
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changes in the parameters. These changes can then be compared with the approxi- 
mate changes obtained by multiplying the entries in Table IV by the appropriate 
changes in the parameters. 


To take a typical example in detail, let us suppose that y, is decreased by 0.1, 
so that Ay,=-0.1. From Table IV and the root values (4.1.2) we infer that 
approximate values for the changed roots will be given by 


q.=0.010 946-0.1 x 0.001 773 =0.010 769 
Gr=4.272 817-0.1 x (-0.022 169)=4.274 986 
r=0.218 118-0.1 x (-0.489 807)=0.267 098 
s=1.415 952-0.1 x (- 0.038 269)= 1.419 779 


5.1) 


Again, using Table III, we calculate the new coefficients of the stability equation as 
A=1 +(-0.1)x(-1)x0 =1 
B=4.72 +(-0.1)x(-1)x1 = 4,82 
C=3.968 +(-0.1)x(-—1)x4.52 =4.420 
D=8.8128 +(- 0.1) x (—1) x 1.464=8.9592 
E=0.096 +(-0.1)x(-—1)x0 =0.096 


and solution of this modified equation gives us the roots 
q.=90.010 772, q-=4.275 081 
r=0.267 074, s=1.418 924 


(5.2) 


Comparison of (5.1) and (5.2) shows no significant error in any of the dampings; 
in the case of the frequency factor s a larger, but still unimportant, discrepancy 
is noticeable. 

The same process has been carried out for each of the parameters involved in 
the equations, with the results shown in Table V. The procedure is identical, except 
for change of i,, and here we have solved exactly both the exact stability equation, 
obtained by proceeding as in (5) of Table I, and the approximate stability equation, 
obtained from (c), Table I. These equations may be put on record: Taking 
Ai, = - 0.005, they are 


O=p* + 4.737 563p* + 4.176 962p? + 8.835 526p+0.096 160: exact 
0=p* +4.73p* +4.17p? +8.8248p +0.096: approximate 


(5.3) 


Table V shows only the exact and approximate root changes due to the changes in 
the parameters specified in the first column. 
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Agreement between the exact and approximate values in Table V is in general 
very good. The damping is given more accurately than the frequency factor, in 
general, and the accuracy of the former is clearly more than sufficient for practical 
purposes, for derivative changes of the orders of magnitude considered, with the 
exception of changes in n,, for which the approximate formule are wildly in error. 
Clearly the method of estimating stability changes suggested here can give excellent 
results. Whether the formule will behave with comparable accuracy for variations 
from another standard set of conditions is open to doubt, but we may 
mention that the case considered was deliberately chosen from the start as one 
likely to be unfavourable. Further comparisons of exact and calculated values will 
be made in a later paper, in which the formule will be used to survey the dependence 
of lateral stability roots upon the minor parameters, over a wide range of 
representative conditions. 


6. ACCURACY OF LATERAL STABILITY ROOTS IN A TYPICAL CASE 


Table V also indicates, for the case considered, the extent to which the 
calculated roots of the secular equation would be altered by changes in the 
individual parameters. If we regard these changes as the probable errors of estimated 
values, then the average probable error in any root may be taken as equal to the 
square root of the sum of the squares of the probable errors due to the separate 
uncertainties. Since these differ considerably in order of magnitude, the overall 
error will not be very different from the largest individual error. 


Clearly, then, the uncertainty in n, is the most important source of inaccuracy, 
affecting particularly the dampings of the important spiral and oscillatory modes. 
Next to it in order of importance comes n,, affecting spiral damping and oscillatory 
frequency: note, however, that the assumed change in n, is larger than the probable 
error of this derivative, which is commonly determined in a model test. An 
important source of error in oscillation damping is uncertainty in y,, which could 
be much reduced if this derivative were more commonly measured in the wind tunnel. 


If these large sources of error could be removed, the stability characteristics 
in this typical case would become much more definite. The assumed changes in 
the other derivatives are larger than the probable errors of estimation, and cause 
no serious errors. In fact, if we could cut down the probable error of y, to 0.03, 
and those of n, and n, to 0.006, it would be possible to specify the roots, with 
their probable errors, as 


»=0.011+0.006, g,=4.3+0.6 
(6.1) 


r=0.22+0.02, s=1.42+0.12 
instead of with errors about four times as great. We may note, also, that the usual 
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neglect of y, and y, is amply justified, unless these derivatives are much larger than 
the assumed changes: on the other hand, i, may have a significant influence on 
oscillation damping. 


APPENDIX A 
Formule for longitudinal stability 


The coefficients of the longitudinal stability quartic are 


m Z mga 1 
Py 


lp 


iz 
x Maw/ at 
ig By lp 


= 1 k | k (—x, tan y+ dt 
Py Py lp lz 


E=k(z.—x, ten»)! — k tan y,)' 


The notation is that of Ref. 7, modified as in the text. We also use the portmanteau 
symbols due to Neumark"), slightly modified : 


N, = 
Kale — Sate 
-(14+ 2) a+ 
R,= —k Xu tan ye) 
Py Py 


— Xw tan Ye) 
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and 


w= — 
(A.3) 
X= — Maws at 
= — p,m, /ig 


The rates of change of the coefficients with respect to the parameters are given 


in Table VI. 


Notation 


APPENDIX B 


The following symbols have the meanings usual in stability theory, as defined 
in Bryant and Gates R. & M. 1801 on Nomenclature for Stability Coefficients. 


in, iz, io 


Fis Be 


S, P; b, m, E, k, 6, 


Xu, Xws Xa, Zus Zws My, Ma, 


Maw/at has been used instead of m,; 


In addition the symbol i, (=4E/mb*) is used, see equation (3.1.2). 


Additional symbols used are as follows : — 


A, B,C, D, E 
Wea Cu Do Bs 
f(p) 


etc. 


> 


coefficients in stability period equation. 
standard values of A, B, C, D, E, respectively. 
stability period equation. 


rates of change of roots with respect to coefficients, the first 
suffix specifying the root, the second the coefficient. Defined 
in equation (2.1.4), formule for calculating and checking in 
Sections 2.2 to 2.5. 


—I,/is, equation (3.1.4). 
ibid. 
— pel, /ix, ibid. 
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N,,P,,Q,, R,S,, T, 
p 
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FLUTTER OF SYSTEMS WITH MANY FREEDOMS 
by 
W. J. DUNCAN, D.Sc., F.R.S., F.R.Ae.S. 


(Professor of Aerodynamics, College of Aeronautics, Cranfield) 


SUMMARY 


Experience has shown that it is often: necessary to retain many degrees of 
freedom in order to calculate critical flutter speeds reliably, but this entails much 
labour. Part I discusses the choice of a minimum set of freedoms and suggests that 
this should be based on the equation of energy and the use of the Lagrangian 
dynamical equation corresponding to any proposed additional freedom. The methods 
for conducting flutter calculations so as to minimise labour are treated in Part II. 


1. INTRODUCTION 


The aim of this paper is to discuss methods for calculating the critical flutter 
speeds and the nature of the motion at these speeds for systems with a large number 
of degrees of freedom. This problem is becoming increasingly important since it 
is now recognised that reliable estimates of critical speeds can, in many instances, 
only be made when many independent kinds of motion of the structure are admitted. 
However, the labour in the calculations increases exceedingly rapidly as the number 
of degrees of freedom is increased. Hence two principal problems arise: — 

(a) The choice of a minimum set of dynamical co-ordinates or degrees of 

freedom which leads to calculated results of adequate accuracy. 


(b) The choice of the method of conducting the calculations after the dynamical 
co-ordinates have been chosen. 


These problems are considered separately in Parts I and II of the paper. 


It is concluded that a particular freedom F must be retained when the balance 
of energy at a critical flutter speed is sensitive to its inclusion, unless it can be shown 
that the amplitude of F is very small. This amplitude will be very small when one 
or both of the following conditions is satisfied : — 


(a) The coupling terms in the Lagrangian dynamical equation corresponding 
to F are all very small. 

(b) The direct impedance for F at the critical flutter speed and for the flutter 
frequency is very large. 


Paper received October 1948. 
{The Aeronautical Quarterly, Vol. 1, May 1949) 
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The equation of energy for any number of freedoms is reduced to a convenient form 
showing explicitly the influence of the phase relations of the motions. It is shown 
that large skew-symmetric components in the aerodynamic stiffnesses will result in 
a large intake of energy when the phase relations of the motions are favourable 
for this. Hence dynamical co-ordinates giving rise to such aerodynamic stiffnesses 
must be scrutinised with special care. 


Some simple illustrative applications given in Section 8 support the validity 
of the proposed method of selecting co-ordinates, but its general usefulness can 
only be judged on the results of extended trials. In any event, it is hoped that the 
analysis of the flutter problem given here will contribute towards a better under- 
standing of the physics of flutter. 


Attention is drawn to the usefulness of comparing the various aerodynamic 
couplings with the corresponding inertial couplings, since this may suffice to show 
that the former are smal! (see Sections 5 and 6). 


Inverse methods appear to be the most advantageous for the calculation of 
critical flutter speeds when there are many freedoms. 


There are a few known special cases where “exact” calculations of critical flutter 
speeds can be made for elastic continuous systems having infinitely many degrees 
of freedom. Such systems throw much light on the general problem of the choice 
of freedoms and they are considered briefly in the following Section. 


2. SOME EXACTLY SOLUBLE PROBLEMS OF THE FLUTTER OF 
ELASTIC WINGS 


Exact solutions of two problems on the flutter of continuous elastic cantilever 
wings have already been published and the details are as follows : — 


(a) Uniform rectangular wing supported by a pair of uniform parallel spars. 
The covering has no proper torsional stiffness and both the flexural and 
torsional stiffnesses are provided by the uniform spars. The inertial and 
aerodynamic coefficients are independent of spanwise position.”: *. *) 

(b) Uniform rectangular wing as for case (a) but now the torsional stiffness 
is provided by the uniform covering or by a torque tube and the flexural 
stiffness is assumed to be negligible* (see the Appendix to Ref. 7). 


Both wings have the following special dynamical characteristics. The modes of 
oscillation in vacuo or in an air stream occur in pairs which, in general, have both 
flexural and torsional components. The displacements in both flexure and torsion 
for the r” pair are proportional to the displacements in the r” mode? of the purely 


* This is an artificial case but it is not without interest since the flexural stiffness does not 
greatly influence the flutter. 


* Or (r—1) overtone. 
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flexural oscillation of a uniform cantilever beam for case (a) and to the displacements 
in the r“ mode of the purely torsional oscillations of a uniform cantilever for case 
(b). When the wings are exposed to an air current the two modes of the r” pair 
become coupled together aerodynamically, but remain uncoupled to all other modes. 
The flutters which occur are therefore effectively binary and the earliest flutter occurs 
in the modal pair corresponding to the fundamental free mode. 


These results are not, in fact, isolated and it is hoped to discuss other similar 
cases in a separate paper. At present we wish to emphasise the remarkable simplicity 
of the solutions, which were obtained in an orthodox way by solving the differential 
equations of the problems. An independent explanation of the results is given in 
Section 6 by means of an analysis of the various couplings between the modes. 


PART I 
CHOICE OF THE DEGREES OF FREEDOM* 


3. TYPES OF GENERALISED CO-ORDINATES 


The generalised co-ordinates should be chosen so as to secure some advantage 
of simplicity or convenience, There appear to be two particularly advantageous 


types : — 
(a) Those for which the specification of the displacements is specially simple. 


Examples. The purely flexural and purely torsional co-ordinates of the 
classical theory of wing flutter. 


(b) Those which secure the greatest simplification of the dynamical equations, 
as by the vanishing of certain coefficients. 


Examples. Normal co-ordinates for which the inertial and elastic coupling 
coefficients vanish. 


In the further discussion we shall adopt normal co-ordinates for the main 
structural deformations. These correspond, strictly, to the free modes of oscillation 
of the structure in vacuo, but differ little from the modes of oscillation in still air. 
They must be supplemented by certain general rigid body freedoms, for the 
resultant force and moment on the system are not, in general, zero when the 
oscillations occur in an air stream. Lastly, there will be one co-ordinate for the 
angular movement of each contro] surface which is concerned in the flutter. The 
main problem is to decide which of the normal modes and general freedoms must 
be retained in the calculations. 


* The ideas of Part I are to some extent based on Ref. 5. 
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For a conventional aircraft, having a fore-and-aft plane of symmetry, and in 
straight symmetric flight, the small oscillations of symmetric and anti-symmetric types 
will be independent. For each type of flutter only the corresponding normal modes, 
general freedoms and control movements must be retained. For example, in 
symmetric flutter a single symmetrically placed rudder takes no part and the body 
freedoms are pitching, normal displacement and fore-and-aft displacement; the 
latter will usually be negligible. 


In order to determine what normal co-ordinates are to be retained the influence 
of an added co-ordinate on the energy balance must be examined. We accordingly 
pass next to the equation of energy. 


4. THE EQUATION OF ENERGY 


For the sake of conciseness and perspicuity we shall use the matrix notation. 
Thus the whole set of the equations of free motion will be written 


Aq+Bq+Cq=0. : ‘ (1) 


where A, B, C are the square matrices of the inertias, damping coefficients and 
stiffnesses, respectively, and q is the column matrix of the generalised co-ordinates. 
A, B and C are independent of the generalised co-ordinates and of time. The 
equation of energy is obtained by pre-multiplying the last by gq’, i.e. by the row 
matrix* whose elements are the generalised velocities, and it is accordingly 


Now A is always symmetric and the kinetic energy of the system is 
a quadratic form in the generalised velocities. Therefore 
dT 
dt 724 49+49 Aq 


since A is independent of time. But 


Aq=4 AG 
by transposition, since A’= A, and we obtain 


* An accent appended to the symbol for a matrix indicates that its transposed is to be taken. 
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For a system subject to aerodynamic actions neither B nor C is, in general, 
symmetric. We may always put 


C=C, +C,, ‘ ‘ ‘ (6) 


where the suffixes i, 2 are applied to symmetric and skew-symmetric matrices 
respectively. Then 


B,=}3(B- B), ‘ (5d) 
and there are similar expressions for C, and C,. Let the dissipation function be 

F=}q' Bg=47' ‘ (7) 


when Y is skew-symmetric and x is an arbitrary column matrix. Also let the 
“potential” be 


so that dt m@C.q. . : (10) 


The equation of energy (2) can now be written 


+ 4 C.q=0. 


The final term in the last equation is of great importance and must be discussed 
in detail. 


The equation of energy is valid for any free motion of the system, but 
particularly illuminating deductions can be drawn when the motion consists of a 
single constituent, corresponding either to a single real root or to a conjugate com- 
plex pair of roots of the determinantal equation. 


Case (a). Real exponential motion (Divergence or subsidence) 


Here we have 
q=qe" ‘ ‘ 
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with A real and qg, a column of real constants. Then 
q = Aq 


by equation (8). Therefore the equation of energy (11) becomes 


If F is negative, obviously the motion can grow. But suppose now that F is positive 
during the motion.* Then we get 


d 


at all instants. So long as V is essentially positive this implies decay of the motion. 
However, V is not wholly of elastic origin and at higher values of the relative wind 
speed it may cease to be essentially positive. Then T may increase, i.e. the motion 
may grow, through dV/dt being negative. At a critical speed for divergence A, F 
and T are zero and the equation of energy becomes 


but the displacements need not be zero as they can be supposed to grow exceedingly 
slowly to finite values. This equation can be made the basis for calculating the 
critical divergence speeds, but we shall not pursue this question here. 


Case (b). Exponential oscillatory motion (Damped, steady or growing oscillation) 


The motion corresponds to the complex pair of roots « +iw of the determinantal 
equation and the displacements can be expressed in the real form 


q =e" (qo. COS wf + q, SiN wf). . (16) 
Hence = + we“ (—q, Sin wt +g, COS wf) 
and C.q=we*™ (—q,’ sin wt + q,’ cos wt)C, (g, COS wt +g, Sin wf) 


= (cos*wt g,'C.g. — sin*wt 


* When the values of the aerodynamic damping coefficients for divergent motion as given by 
two-dimensional vortex sheet theory are used, the condition that the dissipation function 
shall be essentially positive is not satisfied for any value of the non-dimensional divergence 
parameter Ac/V. The dissipation function will be essentially positive when there is enough 
structural damping. 
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for x’Py= —y’Px when P is skew-symmetric, while x and y are column matrices, 
since P’-= —P, Hence the equation of energy (11) becomes 


The special case of this for 3 degrees of freedom is equivalent to equation (42) of 
R. and M. 1155.) We add the explicit expression 


C,=4(C-C’). ‘ ‘ (19) 


Exactly at a critical speed for flutter » is zero and we get by integration of (18) over 
a complete period T 
T 
0 


This is the equation of energy balance at the critical speed. 


5. GENERAL DISCUSSION OF COUPLINGS 


The couplings between the various generalised co-ordinates can be classified 
broadly as follows: — 


(a) Inertial couplings or “products of inertia’ 
It follows from Lagrange’s dynamical equations that the contribution of the 
structure to A,,, the coefficient of g, in the equation corresponding to q,, is 


( Ox Ox dy Oy (21) 
where x, y, z are the rectangular Cartesian co-ordinates of the element of mass ¢m, 
and the summation covers all the elements of mass for which the expression in the 
bracket does not vanish. The complete coefficient A,, contains, in addition, the 
corresponding virtual inertia or reversed acceleration derivative. 

When q, and q, are normal co-ordinates the coefficient A,,=A,, vanishes.* The 
same is true when q, is a normal co-ordinate and gq, is one of the general rigid body 
freedoms. This follows from the fact that motion in any normal mode can persist 
independently of any other motion. 

For a cantilever wing with rigid support at the root we may usually neglect the 
chordwise and spanwise displacements in the flutter. Then the structural contribution 
to A,, is 


* If we assume that the virtual inertias are the same for a finite air speed as for still air then gr 
and qs should be normal co-ordinates for oscillations in still air. 
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Further, when the wing is straight and not swept back or forward and it is assume: 
that each fore-and-aft strip moves as a rigid body,* this expression can be reduce 
to the form (see the list of symbols) 


sArs 


where 


dq, dq; rP dq, dq, + dq, \c; 


(b) _ Cross or compound damping coefficients 


These are mainly or wholly of aerodynamic origin and the array of the 
coefficients is, in general, not symmetric, i.e. B,, ~ B,, in general. We may, if we 
please, analyse these coupling coefficients into symmetric and skew-symmetric parts 
(see equations (Sa), (5b)) and we may call these pure dampings and gyrostatic terms 
respectively. This analysis, however, does not appear to be particularly helpful. 


For the type C wing the general formula for an aerodynamic damping 
coefficient is 


=| 2 Brs (n) dn, (25) 
where 
B. dq, dq, dq. me dq,/ \c, “dq, be 
(26) 


This expression is equivalent to that given in Appendix 3 of R. and M. 1904.“ 


(c) Cross or compound stiffness coefficients 


When normal co-ordinates are used these coefficients are entirely of aerodynamic 
origin and their array is, in general, not symmetric. For a type C wing 


ae 


where, in accordance with R. and M. 1904,‘°) 
dé dt di (c\? 


* To avoid circumlocution we shall henceforward call this a type C wing. 
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6. COMPARISON OF INERTIAL AND AERODYNAMIC COUPLINGS 


When normal co-ordinates are used the inertial couplings are zero and in certain 
instances a comparison of the expressions for the inertial and aerodynamic couplings 
shows that the latter also vanish. More generally, the comparison may enable us 
to show that certain aerodynamic couplings are small. 


The special uniform wings mentioned in Section 2 provide examples of what has 
been said. Here c/c, is unity, the various coefficients are constants, and 


_, do 


29 
da.“ da. 


where k, is constant for each r; this expresses the similarity of the normal and 
torsional displacements in every natural mode. Hence 


di di 


{ mk,k,+ p (kp +k.) +j } (30) 

and all three coefficients vanish when 


7. PRINCIPLES GOVERNING THE CHOICE OF FREEDOMS 


Suppose that a critical flutter speed V. has been calculated on the assumption 
that a particular set of freedoms is effective and let it be enquired whether some 
other freedom, F say, should be included in the calculation. The complete answer 
to the question is given by repeating the whole calculation with the new freedom F 
introduced and by comparison of the new and original critical speeds and modes of 
displacement. This method is usually extremely laborious and it is most desirable 
to find an easier basis of judgment. It is suggested here that judgment may be 
based on a consideration of the possible influence of F upon the energy balance at 
the critical speed V.. If it can be shown that the introduction of F can neither add 
nor subtract more than a relatively very small amount of energy then F may safely 
be neglected. 


The influence of F upon the energy balance will certainly be small when the 
amplitude of F is very small. We may estimate the amplitude of F by forming the 
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additional Lagrangian equation of motion corresponding to F and then assume tha! 
the other generalised co-ordinates have the values already calculated at the critica! 
speed; this approximation will be justified if the amplitude of F is in fact small. The 
amplitude of F will then be given by 


(Sum of coupling terms reversed) / (Direct impedance in F at the flutter frequency). 


Clearly the sum of the coupling terms will certainly be small when all the coupling 
coefficients are small and it may be possible to see that this is so without calculating 
these coefficients. The impedance will be large when the flutter frequency (in the 
absence of freedom F) differs widely from the natural frequency for F alone, 
calculated at the wind speed V , (since the direct stiffness coefficient for F may depend 
on V). We conclude that it will not usually be safe to neglect F when its proper 
frequency is near that of the flutter. On the other hand, when the proper 
frequency of F is much higher than that of the flutter the impedance will be very 
large and the amplitude of F will be very small; hence F can be neglected. 


When the amplitude of F is neither large nor very small the decision to retain 
or reject should be based on a consideration of the sensitiveness of the energy balance 
to F: If some of the stiffness couplings between F and the other degrees of freedom 
show large skew-symmetric components (see Section 4) then F should be retained. 


8. ILLUSTRATIVE EXAMPLES ON THE CHOICE OF FREEDOMS 


A few simple examples are added here to illustrate the principles put forward 
in Section 7. 


(a) In what circumstances can the aileron co-ordinate & be neglected in wing 
flutter? 


The stiffness coupling derivatives L; and H, (zero) are widely different”); like- 
wise M; and Hz are widely different.) Hence the equation of energy shows that 
movement of the aileron may feed in much energy to the system, as is otherwise 
obvious. Consequently the aileron co-ordinates € must be retained unless it can 
be shown that the amplitude in £ is very small. This will be true when either or 
both of the following conditions are satisfied : — 


(i) The aileron impedance is very high. This can only be attained by providing 
the aileron with a very stiff elastic constraint (irreversible control). 


(ii) All the coupling coefficients affecting the aileron hinge moment are 
extremely small. These are the two aileron products of inertia, and the 
aerodynamic derivatives H a,;a:, Haeja:, He. The latter will all be small for a 
certain location of the aileron hinge and the products of inertia can be made 
very small by suitable disposition of masses. 
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(b) Must we bring in overtone modes in considering the flexural-torsional fiutter 
of a cantilever wing not carrying large concentrated masses? 


We have seen in Section 2 that for a uniform cantilever wing with a certain 
special, but quite reasonable, elastic specification, the two gravest modes* are strictly 
uncoupled to the higher modes. For other straight cantilever wings of normal 
proportions the couplings will be small, although not zero. Also, the impedances 
in the overtone modes at the flutter frequency will be large on account of the great 
disparity in frequency. We conclude that the amplitudes in the overtone modes 
will be very small and that only minor errors in the critical flutter speed will result 
from neglect of such modes. 


(c) Flexural-torsional flutter of a cantilever wing carrying a large concentrated mass 

If we take one of the two gravest modes, the inertial coupling coefficient with 
any higher overtone mode is zero. In the expression for this “product of inertia” 
the displacements at the large mass are clearly heavily “weighted,” but these 
displacements are not heavily weighted in the expressions for the aerodynamic 
couplings. Thus we have evidently no right to assume that the aerodynamic 
couplings of the modes are small just because the inertial couplings are zero. Hence 
the higher overtone modes can only be neglected when their impedances at the 
flutter frequency are very large, i.e. when the disparity of frequency is large. Hence 
we must expect to find it necessary to retain some of the lower overtone modes 
beyond the first. 


PART II 


METHODS FOR PREDICTING FLUTTER CHARACTERISTICS 


9. NON-DIMENSIONAL FORM OF THE EQUATIONS OF MOTION 


It is convenient to reduce the equations of motion to a standard non-dimensional 
form. When all the displacements are proportional to exp (Af) the complete set 
of dynamical equations in the selected degrees of freedom can be written in the 
matrix notation 


(AA? + BVA + DV? +E) q=0 (34) 


where A is the inertia matrix 
BV is the damping matrix 
E is the matrix of the elastic stiffnesses 
DV? is the matrix of the aerodynamic stiffnesses. 


* In general both involving both flexure and torsion; the first mainly flexural and the second 
mainly torsional. 
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We shall now suppose that the m dynamical co-ordinates g are all non- 
dimensional (e.g. angles). Then, if p and o are air and material densities and / is a 
typical linear dimension, we may define a set of non-dimensional matrices of 
coefficients by the equations 


where U is one of the velocities of propagation of elastic waves in a typical material 
of construction. Next, introduce the auxiliary non-dimensional variables 


n= JV (p/o) 


The equation (34) can now be rewritten in the non-dimensional form 
(ax? + nbk+d+e)q=0. . (42) 
At a critical flutter speed and for the flutter component we can put 
and the equation to be satisfied is 
(—aw* +inbw +d q=0. . (44) 
The problem consists in finding € and » so that the set of linear scalar equations 


represented by the foregoing matrix equation shall be compatible. Methods for 
doing this are described below. 


10. THE ORTHODOX METHOD OF SOLUTION”) 


We regard » and the matrices a, b, d and e as known constants. If frequency- 
dependent derivatives are used this implies that some estimate of the flutter frequency 
parameter has already been made. 
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The condition for compatibility of the scalar equations represented by (42) is 
the “determinantal equation” 


A (x, £)=90, . (45) 


where the expression on the left is a polynomial in «, in general of degree 2m, whose 
coefficients are functions of €. At a critical value of € this equation is satisfied by 
the pair of equal and opposite roots+iw. Accordingly the parts of A containing 
even and odd powers of » vanish separately and we derive a pair of equations which 
we may write 


A, (w?, &)=0. ; (46) 
The result of eliminating » from these is 


T,-, €)=0. ‘ . (48) 


where T,_, is the penultimate Routhian test function and n is the degree in « of 
the determinantal equation. Equation (48) is a polynomial equation for € and when 
its roots have been found the critical flutter speeds can be obtained from (40). 
Finally, the flutter frequency and the complex modal ratios can be calculated. 


It is worthy of remark that the test function always vanishes when the relative 
density parameter 7 is zero. Thus 7,_, contains a power of 7 as a factor and when 
this is removed we are left with a polynomial in » whose “constant” term is the 
asymptotic form of the test function for low relative air density. Since the critical 
values of € are then independent of 7, equation (40) shows that the “equivalent” or 
indicated critical flutter speeds are independent of the densities. It appears that 
modern heavily loaded aircraft are in or near this asymptotic condition. 


The orthodox method is systematic and yields all the desired information. Its 
drawback is that it becomes unworkably laborious when the number of freedoms 
is greater than three. 


11. THE METHOD OF R. A. FRAZER“) 


Frazer writes the equations (46) and (47) as* 


* Frazer does not use the special reduction to the non-dimensional form adopted here, but 
his method is essentially as stated. 
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where x=-w? ‘ eR) 


Since x and y are both real they could be found from the real intersections of the 
curves (49) and (50). Frazer shows that the full expressions for the functions f 
and g can be calculated from expressions in partial fraction form containing 
the values of A (x, €) corresponding to a special set of points (x, y). These are the 
intersections of a standard set of straight lines satisfying the conditions 


(a) No two are parallel 
(b) No three are concurrent. 


Much labour can be saved by use of this method when the number of freedoms is 
large. It is to be regarded as a special development of the classical method, designed 
to minimise labour. A new lattice has been proposed by Williams?) which is 
thought to reduce the labour still further. 


Frazer and Bratt have proposed mechanical and electrical methods for the 
description of the curves (49) and (50). 


12. METHODS BASED ON THE USE OF MATRICES 


The determination of the roots « for a given € can be made to depend on the 
problem of finding the latent roots of a certain matrix (Duncan and Collar’: '°?). 
An iterative method for finding the latent roots without expanding the characteristic 
determinant is given by Duncan and Collar and other treatments are described by 
Hermann and Dorr."') It is not considered that these methods have yet been 
developed to the stage where they would be advantageous from the aspect of labour 
saving in flutter problems. 


13. THE METHOD OF DUNCAN, COLLAR AND LYON”? 


This method is based on the fact that the mechanical impedance of the system 
is zero at the critical flutter speed and for simple harmonic applied forces having 
the frequency of the flutter. The method consists, in effect, in calculating the 
impedance for a number of assumed air speeds and frequencies and in determination 
of the condition of zero impedance by interpolation. Since the impedance is found 
merely by solving a set of simultaneous linear equations with known coefficients, 
the calculations are relatively simple. The method has the further advantage that the 
modes of oscillation in the flutter are obtained incidentally from the calculations. 
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It is convenient to make the calculations in the following manner. Arbitrarily 
assign the amplitude of one dynamical co-ordinate (say, put g-=1) and use the 
dynamical equations corresponding to free motion in the other co-ordinates to 
calculate the amplitudes of these for the assumed air speed and frequency; this 
requires only the solution of a set of simultaneous linear equations whose coefficients 
are, in general, complex, on acount of the presence of the damping terms. The 
amplitudes of all the dynamical co-ordinates are now known and when these are 
inserted in the equation for the co-ordinate q, the expression on the left hand side is 
equal to the amplitude of the generalised force Q, which must be applied in order 
to maintain the assumed motion. This is the impedance, in general complex, since 
q, is unity. The vanishing of the complex impedance implies two conditions and 
it thus becomes possible to find both the critical flutter speed and frequency by 
interpolation when the impedance has been calculated for a few pairs of values. 


The quick success of this method evidently depends on making a good pre- 
liminary guess or estimate of the critical flutter speed and frequency. 


14. ANOTHER INVERSE METHOD 

We may proceed exactly as in the last method until we reach the final dynamical 
equation corresponding to the assigned co-ordinate g,. However, the expression on 
the left hand side of the equation is equated to zero (free motion) while € in this is 
regarded as unknown. We shall thus derive a value ¢’ for €, in general complex; 
the condition to be satisfied is 


which is equivalent to a pair of real equations. The critical values of € and » are 
to be obtained by trial and interpolation, as in the previously described method. 


LIST OF PRINCIPAL SYMBOLS USED 


The symbols are in alphabetical order and the Greek letters follow the Roman. 
An accent indicates the transposed of a matrix. 

A dot above a symbol indicates its time rate of change. 

The numbered references are to equations. 


A matrix of inertial coefficients. 
Ass general inertial coefficient. See (21). 
a see (24). 
a non-dimensional form of A. See (35). 
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matrix of damping coefficients or, in equations (34) and (36), the 
same divided by V. 


symmetric part of B. See (Sa). 

skew-symmetric part of B. See (5b). 

general damping coefficient. See (25). 

see (26). 

non-dimensional form of B. See (36). 

matrix of stiffness coefficients. 

symmetric part of C. 

skew-symmetric part of C. 

general stiffness coefficient. See (27). 

see (28). 

chord of wing. 

chord at wing root. 

matrix of aerodynamic stiffnesses, divided by V*. See (34). 
non-dimensional form of D. See (37). 

matrix of elastic stiffnesses. 

non-dimensional form of E. See (38). 

dissipation function. See (7). 

rate of change of hinge moment with flexural co-ordinate 9. 
rate of change of hinge moment with rate of change of 9. 
rate of change of hinge moment with torsional co-ordinate 6. 
rate of change of hinge moment with rate of change of 6. 


distance of local flexural centre aft of leading edge as fraction of 
local chord. 


local value of moment of inertia about leading edge, per unit span, 
divided by oc*. 


see (29). 


distance from wing root to reference section or other typical linear 
dimension (y=y//). 


rate of change of aerodynamic flexural moment with aileron angle €. 


rate of change of torsional moment with aileron angle €. 
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l da/dt 


My 


at 


m 


wn wre 


non-dimensional derivative coefficient giving local rate of change of 
lift with normal velocity w. 


non-dimensional derivative coefficient giving local rate of change of 
lift with angle of incidence «. 


non-dimensional derivative coefficient giving local rate of change of 
lift with rate of change of «. 


non-dimensional derivative coefficient giving local rate of change of 
moment about leading edge with normal velocity w. 


non-dimensional derivative coefficient giving local rate of change of 
moment about leading edge with angle of incidence «. 


non-dimensional derivative coefficient giving local rate of change of 
moment about leading edge with rate of change of <. 


number of dynamical co-ordinates, in Sections 9 and 10. 
local value of mass per unit span, divided by oc’, in Section 5. 


local value of mass moment about leading edge, per unit span, 
divided by oc’. 


column matrix of dynamical co-ordinates, which are taken to be 
non-dimensional. 


see (16). 

r” dynamical co-ordinate or its amplitude. 
kinetic energy of system. 

periodic time of flutter. 

penultimate test function in Routh’s sequence. 


time. 


a velocity of propagation of elastic waves in a typical material of 
construction 


“potential” of system. See (9). 

true air speed. 

lateral distance from plane of symmetry or wing root. 
normal displacement. 

y/l in (23), (25) and (27). 

v‘(p/o), in Sections 9 and 10. 

local twisting or pitching displacement. 

downward normal displacement of leading edge, divided by c,. 


aileron angle in Section 8. 
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€ non-dimensional speed coefficient elsewhere. See (40). 

K a non-dimensional frequency parameter. See (41). 

d coefficient of ¢ in exp (AP). 

m real part of A. 

p air density. 

o density of typical material of wing. 

w coefficient of imaginary unit in A or x. See (16),-(43). 
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NOTE ON PROPELLER-TURBINE 
REDUCTION METHODS 
by 
E. C. PIKE, M.A., A-F.R.Ae.S. 


(Member of the Technical Staff of the Royal Aeronautical Society) 


1. INTRODUCTION 


It has been shown") that since the performance of a turbo-jet engine 
can be expressed theoretically as a function of the two variables N/ /@ and V/ V4, 
the performance of a turbo-jet aircraft can also be expressed as a function of two 
variables, N/ /6 and W/&. Two methods can then be used for the reduction of the 
performance to standard conditions. In the “non-dimensional” graphical method, 
flight test measurements of, for instance, V/ / 6 are plotted against N/ /@ and W/6, 
and the standard value of V/¥/@ is read from these curves at the standard values 
of N/ /@and W/6. In this way, the relation between the parameters is established 
during the tests. The alternative method is an analytical one in which corrections 
are applied to, for example, the level speed for differences from standard in the 
independent variables. For this method the effect on the performance of changes 
in these variables must be predetermined. 


The case of the propeller-turbine aircraft is more difficult than for turbo-jets 
for two reasons : — 

(i) A third variable is introduced into the engine performance equations, which 
makes a “non-dimensional” graphical method awkward from the flight-testing view- 
point. 

(ii) There is little information on the.variation of engine performance with the 
independent variables, because little test flying has been done, and there is no 
variable altitude test plant in this country. 


For propeller-turbine aircraft there are three possible ways of reducing the 
performance to standard conditions : — 

(a) Bya “non-dimensional” graphical method (with an extra variable), in which 
the relations between the variables are established during the test, as for turbo-jets. 

(b) By a derivative method, in which the variation, for instance of F/6 with 
V/ J, is established previously by test. 

(c) By a derivative method, in which the variation, for instance of F/& with 
V//9, is assumed from firms’ estimates. 


The possibilities of establishing such methods are discussed. 


Paper received June 1948 
{The Aeronautica) Quarterly, Vol. 1, May 1949) 
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2. NOTATION 


Temperature and pressure symbols without suffix denote atmospheric conditions 
Other suffixes to temperature and pressure symbols refer to the positions indicated 


in Figs. 1 and 2. 


4 
9 
W 
a @ 
= 
3 


Fig. 1. 
Single shaft engine. 


f. (V | V6) 
brake horsepower supplied to propellers (h.p.). 
f. ow 
propeller power coefficient. 
aircraft total drag (Ib.). 
Of, 
f, 
net jet thrust (Ib.). 
total propulsive thrust power (ft. Ib. / sec.). 
0 (BHP 4) 


BHP/iV6 

F/é 
propeller advance ratio. 
V/V0 On 


mass flow (lb./sec.). 
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JET PIPE 


| COMPRESSOR 
\ TURBINE 
TURBINE 


| 

| 


Fig. 2. 
Free turbine engine. 


N engine speed (r.p.m.). 
No compressor rotational speed (r.p.m.) 

free turbine engines. 
Np propeller shaft rotational speed (r.p.m.) } 


air pressure (Ib. /in.*). 
Q fuel flow (Ib./sec.). 
T absolute temperature (°K). 


THP total thrust horsepower. 

V forward speed (knots). 

equivalent air speed =V (knots). 

Te equivalent air speed for minimum drag (knots). 

W aircraft weight (Ib.). 

= (as in Ref. 1). 

A’(W/6) correction to W/®=(W/®),- | and similarly for 
A(W/s) AV, A’V ete. 
propeller efficiency. 

6 relative temperature = T /288. 

VilVing° 

& relative pressure = p/ 14.7. 

relative density. 

MY) the “third variable” (see Section 3). 
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3. THE THIRD VARIABLE 


An ordinary turbo-jet unit works on a single operating line such as AA’ (Fig. 3) 
(Actually this line fans out for different values of V // @ below A’, where condition 
are not choking.) This means that for any value of N/ 4 the fuel flow is defined, 
and the “non-dimensional” performance parameters can be expressed as function 
of two variables only, (N//6, V/ 


In the case of a propeller-turbine the operating point can move round the 
“constant N//@ line” (BB’ say), because the distribution of power between the 
propeller and jet thrust can be varied. This corresponds to a variation in T,/T or 
propeller pitch angle, or whatever is chosen as the third variable. 


Broadly speaking, propeller-turbines can be divided into two classes, “single 
shaft” engines with the propeller geared to the compressor-turbine shaft, and free 
turbine engines in which a separate power turbine drives the propeller. S. J. Moyes 
of The National Gas Turbine Establishment, in an unpublished note has shown that, 
for constant exit area, in the first case, 


BHP/(é¥v and Q/(®¥ 8) are all functions of N/ /0, V//6 and T,/T 
and in the second case 
BHP/(®/) and Q/(®¥ 8) are all functions of N,/ / 0, V//@ and N,/ 


In either case the third variable could equally well be pitch angle or T,/T, but 
it should be noted that in the second case (free-turbines), N,//@ has replaced 
T,/T because it is possible to have a variation in Np//@ at constant N,//@ 
and T,/T. 


It has been found that firms prefer to use jet pipe temperature 7, instead of 
turbine inlet temperature 7, in the third variable, because it is easier to measure. 


If the exit area is not fixed, a fourth variable is introduced, unless the area can 
be expressed as a dimensionally correct function of the existing variables. For 
example, it might be arranged that the exit area was a function of Mach number 
G.e. of V//6). The case of a fourth independent variable will not be considered 
here, as it does not occur among current engines. 


In some cases the three fundamental variables may be reduced to two by a 
suitable control arrangement which must be dimensionally correct. For example, in 
the case of a single shaft engine if 7,/T or pitch angle is related to N/ 6 and 
V/ 4/0, there are only two independent variables, but if, say 7, is kept constant, 
there will still be three variables. Again, in the case of the Bristol Theseus, N../N> 
is kept constant (but not by a direct gearing), so there are only two independent 
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Fig. 3. 
Compressor characteristics. 


variables (for constant exit area). For a free turbine engine, a simplificatian arises 
for fuel consumption when the power turbine is choking. In this case 


9) =f (No/ V9, V/V 9) 
irrespective of any relation between N, and N,. 


The third variable might be eliminated under certain conditions—for instance 
the propeller pitch angle might be constant in the cruising range—but the possibility 
is remote. 


4. OTHER VARIABLES 


Combustion efficiency may vary with pressure; calorific value of the fuel may 
alter; and there may be viscosity effects on the engine performance. None of these 
variables is included in current turbo-jet reduction methods, but it has been 
suggested that viscosity effects may be important. 
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5. PROPELLER EFFICIENCY 
It is shown in Appendix I that 
n=f, (N/V, ¥). 
where vy is the third variable, (7,/T, T,/T, or Np/ V9). 
Therefore, if a graphical method is used, the propeller efficiency will not enter 
the problem as a separate item, as shown below. If an analytical method is used, 


the variations of » with the independent variables must be determined separately 
and will depend on the variations of BHP /(#1/ ) with these variables. 


6. POSSIBILITIES OF A “NON-DIMENSIONAL” GRAPHICAL METHOD 


Assuming the engine performance to be a function of three variables, we have: 


F/5=f, (N/V9, V/V9, wv) . (1) 


(where ¥ is the “third variable”), and it is possible to use a graphical method, 
although owing to the presence of the third variable, the flight test programme 
required is considerable. 


For example for level speed : — 
THP = DV /326 


THP DV 


l 
326 
BHP FV 1 DV | 


i.e. 


or 326 V9, . (4) 


since D/i=f, (V/V 6, W/®) as for turbo-jets. Therefore in level flight 


V/V6=f, (N/ V9, 
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Therefore if curves of V//@ against N//@ are plotted for several values of 
W/s and y¥ under ambient conditions, standard values of V//@ can be 
read from them. 


Climb (assuming a relation between V//@ and as for turbo-jets‘)) 
and range can be dealt with in a similar way and are functions of three variables, 
except in the special cases where the engine performance is a function of V/ 4 
and one other variable only, when the method is the same as for turbo-jets. 


7. POSSIBILITIES OF AN ANALYTICAL METHOD 


The difficulty here is that the relations between the engine performance and the 
variables are not established, and that in addition the variation of propeller efficiency 
must be included. Assuming that these relations can be established, there are 
several bases on which the performance can be reduced. For example, for level 


speed : — 
THP = DV /326 


326 THP/ 0)=(D/8) (V/V 9) 
and as in (4), the left hand side is a function of (N//6, V/@, wv). Therefore 
(V//9)(D/o)=f, (N/V 9,V V9, (5) 


where f, represents the total “non-dimensional” thrust horse power, i.e. 


(6) 
It is now assumed that 
p=avs 


and therefore, as shown in Appendix II 


a(D/a) 


W/® 2 
D/i a(W/s) 
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7.1. Differentiating equation (5) at constant N/ 8 


of, 


fi fe 
V//0 of. 
where B= i.e. the percentage change in thrust horse- 
ee power for one per cent. change in V/ 6 or 
of. in 
f, ov 


Using the fact that (V/ /@) (D/)=f, and dividing through by f,, 


6) ct 


d(V/¥/8) Z\d(W/a) 
(7-5) W/o Vivo 


V/J/0 


Rearranging, and writing as a correction 


A’ V/V 9) Z\A(W/5) 


Therefore it would be possible to correct the results at constant N/ ¥ 4 which 
is the most convenient method for turbo-jets, but the difficulty arises in the deter- 
mination of B and C. (It should be noted that reduction at constant N/ ¥ 6 neglects 
change in drag due to Mach number.) 


To find the form of B and C:— 


Differentiating equation (6) partially with respect to V/ v6 


a( 
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where K = 
n O(V/¥/8) 
. i.e. the percentage change in », BHP or 
G ive eee thrust (“non-dimensional”) for one per 


~  d(V/ 


H= 


cent. change in V/ 6. 


Simplifying, 
B=K+G+(1+H-G-K)x(F/8) (V/V 9)/f, 
C can be expressed in a similar way in terms of the variation of », BHP /(é¥ @) and 


(F with v. 


7.2. Atconstant V//é 
A(N/ V4) Av 


here E= N/V 9 __ 
where E f, 


and is about as complicated as B. 


The advantage here would be that change of drag due to Mach number is 
avoided in using constant V/  @. 


7.3. At constant temperature 


W @ 


A similar expression can be derived at constant pressure. 


8. DATA REQUIRED TO USE AN ANALYTICAL METHOD 


It is obvious from the foregoing that in order to apply an analytical method, it 
is necessary to know the percentage change in 


n, BHP /(6/ 9), and Q/(6/ 4) 
due to one per cent. change in any two of the three variables. 


N/ V8, ¥. 
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The change in » with the other variables could possibly be determined in a 
general way from available data on propellers, and it would probably be necessary 
to use different values for cruising, top speed and climb conditions. 


The changes in the other parameters could be left for the individual operator 
to fit in the values applying to the particular aircraft-engine combination (e.g. 


N/V9 OF 
X= - jets) 
F d(N/¥V8) for turbo-jets®)). It would be more satisfactory if 


representative values could be determined either from flight tests, or from firms’ 


estimates, but it seems unlikely that this is possible. 


9. CONCLUSIONS 


A graphical method of reduction is much simpler than an analytical method, 
provided enough flying can be carried out to cover three variables. In special cases 
where there are only two variables the graphical method is always easier. 


An analytical method is feasible, but hampered by lack of generalised informa- 
tion on the variation of n, BHP/(é/ 4), F/® and Q/(6¥ 4) with the variables N/ ¥ 4, 
V//6,w~. Reduction at constant N//@ or constant V//@ seems simpler than at 
constant temperature or pressure. 
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APPENDIX 1 
n=f (Cp, J, Mach number). 


BH 
But C, varies as ashe for constant diameter 


(N/V Op since 
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V/J/0 
J varies as N/V for constant diameter 


Mach number varies as V/ / 6 


and therefore, since =f,(N/ V9, wv) 


n=f,(N/ V9, V/ V9, 


APPENDIX Il 


2 


, bW 
D=aV,? + —— where a, b are constants. 


rr 


Therefore 


| (YL) 


Differentiating partially with respect to V/ /é 


=2aV 


| vive 


METHODS 


6 


Djs | 


aV? - bW?/V? 
aV? + bW?/V? 


D_ 
=Z 
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Differentiating partially with respect to W/é 


a(D/s) W 
[vive 
aV?+bW?/V? 
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DETERMINATION OF DRAG OF 
JET-PROPELLED AIRCRAFT FLIGHT 


by 


G. W. TREVELYAN, M.A., 
and 
D. R. BLUNDELL, Grad.R.Ae:.S. 


(Aerodynamics Department, de Havilland Aircraft Co. Ltd.) 


1. INTRODUCTION 


In order to obtain full knowledge of the performance of an aircraft, one of the 
important quantities to be determined is the drag. The usual procedure is to measure 
the level speed of the aeroplane and to calculate the drag under these conditions 
from the engine thrust, as obtained from the engine makers’ specification. To do 
this, accurate thrust data are essential and in certain cases these may not be 
readily available. 


In order to eliminate the possibility of systematic errors in the determination of 
engine thrust, a method has been developed whereby the drag at low Mach number 
can be measured directly under such conditions that accurate knowledge of the engine 
performance is of minor importance. This method is satisfactory only on jet- 
propelled aircraft, since on conventional aircraft there are the added complications 
of propellers, slipstream, radiators and so on. 


The method has been used to determine the drag of the de Havilland D.H. 100 
Vampire and has given good results. The drag determined from earlier tests by 
the usual procedure is found to be in good agreement with that given by the new 
method. 


2. THE METHOD 


The motion of an aeroplane in straight flight is defined by the following equation 
based on Newton’s Second Law (see Notation): — 


at 


When the aircraft is flying with insufficient thrust for sustained steady level spzed, 


F=C)Sq+ (1) 


Paper received October 1948. 
[The Aeronautical Quarterly, Vol. 1, May 1949) 
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HORIZONTAL 


Ww 


Fig. 1. 
Force diagram. 


straight flight can be maintained either by decelerating at constant height or by 
gliding (at constant or varying forward speed). 


The method consists of either decelerating or gliding at minimum r.p.m. and 
determining the deceleration or angle of glide. By this procedure, the thrust term 
is a relatively unimportant factor in equation (1), and the drag deduced is not 
dependent on accurate knowledge of the engine performance. 


In deceleration at constant height, sin 6 is zero, and the important quantity to 


be measured is the deceleration, -dV,/dt. It is assumed that, over the distance 
covered, the wind velocity is constant, so that 


dV, _ dV; 
dt dt 


(2) 


Modifying equation (1) and including equation (2), the following expression 
is obtained for extra-to-induced drag at 100 ft./sec. at sea level: - 


dhs) 


(3) 


In a glide, on the other hand, the true speed is ideally kept constant, the 
deceleration therefore is zero, and sin 6, which is obtained from the rate of descent 
and forward speed, becomes the important quantity. In practice it is easier to 
maintain constant indicated speed, and consequently the analysis becomes more 
complicated. For this reason, and others to be discussed later, more emphasis has 
been placed on the deceleration method. 
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w cos > 


HORIZONTAL 


Fig. 2. 
Triangle of velocities. 


2.1. Theory 


Figure 1 shows diagrammatically the external forces on an aeroplane in straight 
flight in a vertical plane. It is assumed that the thrust line is parallel to the 
relative air flow. 


Figure 2 shows the corresponding triangle of velocities, taking into account 
the wind velocity. 


The diagrams show the aeroplane in climbing attitude to facilitate the use of the 
usual sign convention. 


The derivation of equation (1), showing the approximations made therein, is 
given below: — 


Resolving parallel to the relative air flow, by Newton’s Second Law, 
Wd 


It is now assumed that the cosines of the angles 6, ¥ and « are unity: so 


; : : (5) 
and 
L=W=C,S¢q . ‘ (6) 
Using equation (6), the induced drag is then given by: — 


wi 
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Also, by definition, 
and 
Equation (4) then resolves into the form of equation (1). 


Now, referring to Fig. 2, 


dh 1 
sin 6 = dt V; ° (10) 


The equation of calibration of the altimeter, based on hydrostatic laws, is 
as follows : 


where 5h, is the change in altimeter reading for a change of pressure, dp. 


The change, 5h, in true or “tape line” height for a change of pressure dp 
is, however, given by: — 


From equations (11) and (12), 
dh _o dh, 
2 
os 
= 


Substituting for dh/dt in equation (10), 


In the decelerations, 6=«=0. This method therefore is potentially more 
accurate, since half the approximations of the generalised theory are eliminated. 
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If, for some reason or another, constant altitude is not maintained, a correction 
of the following form can be applied : — 


AD . Wsin 6 
q 


3. DESCRIPTION OF THE TESTS 


Tests were made on a D.H.100 Vampire, fitted with a Goblin 2 engine. The 
aircraft was in standard service condition, except that the gun ports were faired over. 
No special cleaning or polishing was done. 


The aircraft was fitted with an automatic observer containing altimeter, air 
speed indicator, ambient air thermometer, time base, engine r.p.m. indicator, jet-pipe 
total-head gauge and jet-pipe temperature gauge. These were photographed by a 
camera operated as required by the pilot from a button on the control column. 
Ambient air temperature was also obtained from the Meteorological Office, 
Dunstable, as confirmation of that from the electrical impact type thermometer 
transmitting to the automatic observer. Wind data were also obtained from the 
Meteorological Office. 


As stated before, greater emphasis was placed on the level flight decelerations 
than on the glides because of the easier flight programme for the pilot, the simpler 
reduction and greater accuracy of the results and also, the wider speed range that 
could be covered. 


For the decelerations, the aircraft was flown level at 10,000 ft. at approximately 
the Mach number at which the compressibility drag rise sets in. The throttle was 
closed rapidly, reducing the r.p.m. to the minimum, and the aircraft held in straight 
and level flight until the speed had dropped to 140-150 m.p.h. 1.A.S., camera shots 
being taken at frequent intervals. 


Allowance was made in the flight programme for the effect of position error 
on the altimeter, and the accuracy of the piloting was such that the true pressure 
altitude was maintained in each case within +80 ft. from a mean value. 


The glides were made at minimum engine r.p.m. at constant LA.S. from 
10,000 ft. to 1,000 ft. 


Four decelerations and two glides were done, the whole of the tests consisting 
of two flights of half an hour each. 


Tests were also made on other aircraft with equally satisfactory results. 
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TRUE AIR SPEEDO ~ M.P.H 


° 20 40 60 80 100 120 140 160 180 200 
TIME FROM START OF RUN ~ SECONDS 
Fig. 3. 
True air speed and time from deceleration No. 4. 


4. ANALYSIS OF THE TEST RESULTS 
The deceleration tests were analysed as follows : — 


The weight at any time during the tests was determined by estimates of fuel 
used and flying time and was checked by fuel gauge readings and by the total fuel 
used during the flight. 


From each shot of the automatic observer film records, the true air speed was 
determined, and these values were plotted against time throughout the test. Fig. 3 
gives an example (Deceleration No. 4) of the True Air Speed—Time curve obtained, 
showing the mean line used in the analysis. 


The deceleration, —dV,/dt, and hence the inertia term, —(W / (dV ,/dt), were 
obtained by graphical differentiation of this mean curve. Although graphical 
methods are not in general fully satisfactory, they do appear to be the best available 
for this purpose. Since the curve at the extreme beginning and end of the decelera- 
tion cannot be determined accurately, two or three shots at either end have been 
omitted in each case. 
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The net engine thrust (including intake drag) was obtained, as usual, by 
determination of engine mass flow and jet velocity from the readings of jet-pipe 
total head and temperature. 


The induced drag, D,, was evaluated on the assumption that the factor, k, 
was 1.15 (see equation (7) and Section 5). 


Knowing these quantities, D,,,, the extra-to-induced drag at 100 ft./sec. at 
sea level, was obtained from equation (3). 


Table I shows the main figures of the calculation for three representative shots 
from the same test as Fig. 3 (Deceleration No. 4). 


TABLE I 
Time from start t sec. | 17 68 143 
True press altitude ft. | 10,060} 10,040 | 10,000 | 
| | | 
Ambient air temperature ¢. +7 | +7 | +7 | 
True air speed m.p.h. 416 | 300 | 198 
| Dynamic pressure q Ib. / sq. ft. 312 163 | 71 | 
| Weight W ib. | 7,770) 7,770| 7,770 | 
Deceleration - ft./ sec.” | 4.49 2.48 | 1.79 | 
| | 
Inertia term ie Ib 1,082 | 598 | 431 | 
@ | 
Thrust F Ib. | 46 35 | 31 | 
Induced drag Ib. 44 85 195 | 
Extra-to-induced drag CySq Ib. 1,084 548 267 | 
Extra-to-induced drag at | 
100 ft./sec. at sea level Pixs Ib. | 41.4 40.0 44.8 
Reynolds number - R | 18.7 x 10° |13.5x 10° | 8.9 x 10° 
= 
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The drag from the glides was obtained in substantially the same way, except that, 
in addition, graphical differentiation of the True Pressure Altitude—Time curve was 
necessary to give R, and hence W sin 6 (see equation (14)). Also in order to take 
account of any variation of wind speed with height, 


d 
cos 4) 


was substituted for —- dV ,/dt in equation (3). 


The values of D,,, deduced from all the decelerations and glides are plotted 
against Reynolds number in Fig. 4. 


Below a Reynolds number of about 16 x 10°, D,,,, shows the expected behaviour, 
increasing with decreasing R. Previous tests on the Vampire (see Section 6) 
have indicated that D,,, remains sensibly constant over a wide Reynolds 
number range from 12 x 10° to 27x 10°. From these deceleration tests, the value 
of this constant drag is about 41 Ib. 


It would appear, therefore, that the glides (particularly Glide No. 1 at 
R=24x 10°) do not give as satisfactory results as the decelerations, and the latter 
method is therefore recommended as the better and more accurate. 


5. ACCURACY OF THE METHOD 


The aim of this method is the determination of the drag (at low Mach number) 
directly from flight tests under conditions when any systematic errors in engine 
performance have only a small effect on the result. 


It will be seen from Table I that this aim has been achieved, since an error 
of +10 per cent. in the thrust causes not more than about + one per cent. alteration 
to the drag. With the usual procedure of level speed tests, this error in the thrust 
would cause slightly more than 10 per cent. change in the drag. 


The deceleration, —dV/dt, as determined by graphical differentiation, is con- 
sidered to be accurate within +5 per cent. A test was made using a sensitive 
accelerometer (by the method developed at the Royal Aircraft Establishment for 
high Mach number dives), but it was found that, due to the large change of incidence 
during the run, the corrections for the attitude of the instrument were excessive and 
the results unsatisfactory. 


Allowance cannot be made in this analysis for the effect of horizontal and 
vertical air currents, since it is impracticable to collect all the requisite meteorological 
data. An up-current of 2 ft./sec. would cause a decrease in the deduced drag of 
about 2 per cent. at the highest speed and about 20 per cent. at the lowest. The 
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tests were, however, all done on calm, fine days, and each pair of decelerations was 
made over a different district, covering in each direction about 15 miles in about 
three minutes. By taking these precautions, it is thought that the effects of an: 
gusts encountered were eliminated. 


Examination of Fig. 4 shows a scatter of about +6 per cent. in the test 
points. However, if a statistical analysis is made of the results available from the 
decelerations (97 test points) taking into account all the random errors in 
deceleration, weight, speed and so on, it appears that the accuracy of the mean line 
is about +1.5 per cent. (95 per cent. probability). 


Beyond these random errors, there may be a systematic error due to the use 
of an incorrect value of the induced drag factor, k. However, from analysis of 
these and other tests on the Vampire, the true value might be expected to lie within 
+5 per cent. of the assumed value, kK=1.15. This error would have negligible 
effect on the deduced drag at the higher Reynolds numbers, but would cause an 
alteration of about +4 per cent. at the lower speeds. 


6. COMPARISON WITH OTHER TESTS 


In level speed tests on a number of other Vampires, the thrusts in flight were 
determined from the readings of jet-pipe total head and temperature and, using 
these, the drag was deduced in the usual way. Over the Reynolds number range 
of these tests (11 x 10° to 28 x 10°), it appeared that D,,, for the Vampire was 
constant at about 42 Ib., showing good agreement with the results from the 
decelerations. 


Level speed tests were also made with and without transition wires at’ 4 per 
cent. chord on both surfaces. These tests showed that the wires caused a loss 
of speed of 23 m.p.h. at R=12x 10°, decreasing to five m.p.h. at R=27x 10°. 
This indicates that the effective transition point moved forward from about 20 per 
cent. chord at the lower Reynolds number to the leading edge at the higher. This 
change in the effective length of laminar flow is sufficient to explain the fact that 
Do. Temains constant over this Reynolds number range. 


7. CONCLUSIONS 


The drag (at low Mach numbers) has been determined directly from flight tests 
under conditions when the engine thrust is an unimportant factor. 


Results from a Vampire show, in the Reynolds number range from 21 x 10° 
to 16x 10°, that D,,, is constant at about 41 Ib. at 100 ft./sec. at sea level. Below 
this range the drag increases to about 45 Ib. at R=8x 10°. Although there is some 
scatter in the test results, the accuracy of the mean curve is considered to 
be +1.5 per cent. 
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It appears that the deceleration method is better and more accurate than the 
glide method. 


Tests have been made on other aircraft, with equally satisfactory results. 


The drag obtained by this method is with small mass flow through the engine. 
With greater mass flow there may be differences in the external flow over the intakes, 
and allowances may have to be made in applying these results to conditions 
of high thrust. 
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NOTATION 
F total propulsive thrust. 
L total lift. 
W instantaneous aircraft weight. 
D total aircraft drag, = D,+D,. 
dD, induced drag = : ‘ 
D; total extra-to-induced drag =C)Sq. 
Beas extra-to-induced drag at 100 ft./sec. at sea level. 
k induced drag factor. 


wing span. 


gross wing area. 


Ci coefficient of lift, based on gross wing area. 

Cp coefficient of extra-to-induced drag, based on gross wing area. 

6 angle between relative air flow and horizontal (positive for climb). 

v angle of drift in horizontal plane (usual definition). 

€ angle of drift in vertical plane, as defined by Fig. 2. 

9 angle between wind direction and aircraft headings (=wind direction 
— true course). 

V, velocity relative to ground. 

Vy true air speed. 

V; equivalent air speed. 


b | 
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indicated air speed. 
horizontal wind speed. 
true altitude. 


pressure altitude. 
true or “tape-line” rate of climb. 


dh od 
= =rate of change of pressure altitude (positive for climb). 


Reynolds number. 
time. 


dynamic pressure=4p,V;* (suffix ,,, refers to 100 ft./sec. at sea level). 


Atmospheric quantities 


T 


atmospheric density. 
relative atmospheric density. 
atmospheric pressure. 


atmospheric absolute temperature. 


Suffix , refers to sea level conditions. 
Suffix , refers to ambient conditions at pressure altitude. 


Suffix , refers to I.C.A.N. conditions at pressure altitude. 
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NOTES ON THE LINEAR THEORY OF 
INCOMPRESSIBLE FLOW ROUND 
SYMMETRICAL SWEPT-BACK WINGS 
AT ZERO LIFT 


by 


F. URSELL 


(Department of Mathematics, The University, Manchester) 


SUMMARY 

In much of the recent work on the flow of an incompressible fluid past a 
swept-back symmetrical wing of small thickness at zero incidence it has been 
assumed that the disturbance introduced by the wing can be represented by means 
of a distribution of sources and sinks over the centre plane of the wing, the source 
strength being determined from the equation of the wing profile. It was found by 
Neumark") that on the wing itself this procedure was liable to lead to serious 
errors. Neumark and Kiichemann”) studied certain special types of wings 
for which they were able to find the correct value of the pressure at the surface by 
evaluation of the pressure at points off the wing and continuity considerations. This 
process is rather cumbersome; it may therefore be convenient to regard the problem 
from a rather different point of view. 

In this paper a derivation is first given of the linearised equations of motion. 
Green’s theorem on potential functions then leads at once from these equations to 
the method of sources, whereby the pressure is expressed as a surface integral 
extended over the wing. For points on the wing, the principal value in some sense 
of the surface integral must be taken. To this must be added a line integral 
depending on the definition of the principal value. When the swept-back wing 
consists of two equal cylindrical surfaces joined at an angle, an integral equation 
can be derived giving the pressure distribution along the centre section in terms of 
the profile, agreeing with that given by Neumark and Kiichemann. 


Paper received July 1948. 
(The Aeronautical Quarterly, Vol. 1, May 1949) 
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It is next shown that this equation can be solved explicitly) without 
making use of the numerical process devised by Kiichemann. The solution is 
applied to the simple pressure distribution studied by the latter, and it is confirmed 
that the profile corresponding to this distribution crosses itself. A closer 
examination of the solution shows that this difficulty can be avoided, if the pressure 
distribution is slightly modified near the edges, where, on the linear theory here 
considered, it becomes infinite. (It is known that the linear theory is incorrect near 
the edges, so these infinities are not disturbing.) 


When two cylindrical surfaces are joined at an angle to form a swept-back 
wing, the wing plan, consisting of straight lines, has a discontinuity in gradient near 
the nose. Suppose that the wing plan is slightly changed so as to have a small but 
finite radius of curvature at the nose, the section remaining unchanged. It is shown 
that the change in the pressure distribution along the centre section tends to zero 
as the radius of curvature decreases to zero, so that the discontinuity in the gradient 
of the wing plan has no pronounced effect on the pressure distribution. 


It is next shown that it is not possible to design wings with kinked isobars in 
the centre, unless suction slots or similar devices are used. Finally, it is proved 
that there may be no wing or one wing with a given plan form and pressure 
distribution, but that there cannot be more than one. 


1. DERIVATION OF LINEAR THEORY 


Consider the flow past a wing possessing a central plane of symmetry, which 
is taken as the plane y=0 (Fig. 1); the plane z=0 is taken perpendicularly to y=0 
through the chord of the centre section; the plane x=0 is taken bisecting the chord 
of the centre section, the nose being at x= +1. 


Suppose that the wing is symmetrical about z=0 and at zero incidence to the 
stream. Then no trailing vortex sheets are formed in the fluid and a velocity 
potential ® exists, such that the velocity vector grad ® is continuous in any region 
excluding the wing. 


Write 
P=-V.x+9. ‘ ‘ (1) 
where V, is the velocity of the main stream; ¢ satisfies 
do 
and grad ¢->0 as the distance from the wing tends to infinity. . , ‘ (3) 


The whole motion is symmetrical about z=0; it is therefore sufficient to consider 
only the region z > 0 in this paper. 


Take as the equation of the wing 
z=(G,))> 0; .. ‘ (4) 
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The region of the plane z=0 where ( (x, y) > 0 is called the wing plan A. Outside 
A, satisfies 4¢/0z=0 on z=0. There is no velocity component normal to the 
wing, i.e. 


It is assumed that if ¢ (x, y), is of the form eZ (x, y), where ¢ is a small parameter, 


@ may be expanded in powers of «. A linear theory is obtained by ignoring powers 
of e higher than the first. 


+! 


Fig. 1. 
Co-ordinate system. 
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The linear potential satisfies 


of =0 when z=0 and (x, y) is on A. (5) 


ox + oz > 


~ 


0@/0z is thus prescribed on the whole plane z=0 when ( (x, y) is given. Also 
vanishes at the leading and trailing edges, whence 


L.E. oo L.E. yt 
| a dx=-V, | - =0 when z=0 and y has any constant value. (6) 
T.E 


The pressure on the wing is 


constant 4p | grad | =constant+ pV, to the first order. . (7) 


The components of the vector grad @ will be denoted by (u, v, w). 


2. CALCULATION OF THE PRESSURE ON THE PROFILE 


The linear equations (1) to (7) form the basis of the linear theory. It has been 
seen that » is symmetrical about z=0, i.e. 


Oz 


~ 


except on the wing plan A where 


Od Xx 


and is a given function. 


If the wing is finite, grad » vanishes for large values of x*+y*+2z", and the 
normal derivative is given on z=0. With these conditions » is uniquely defined in 
z> 0 as the potential due to a continuous distribution of sources over the wing 
plan; the pressure is proportional to the corresponding x-component of the velocity. 
The double integral defining ¢ becomes meaningless when z is put equal to zero, 
and near the plane z=0 the convergence is non-uniform. This was overlooked by 
Ludwieg”), whose calculation was corrected by Neumark. Neumark used 
the continuity of 0¢/x with respect to z: 


op _ lim (**) 
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The latter quantity is uniquely defined and can be evaluated for the class of cylin- 
drical wings considered by Neumark. 


It will now be shown from Green’s formula that the definition of u as a double 
integral still holds on A, provided that a line integral depending on the definition 
of the surface integral is added. 


2.1. Application of Green's Theorem 

It is thus required to determine u(X, Y) in terms of w(X, Y) on the wing. 
Denote by R the position vector of (X, Y, 0), and by r the position vector of 
(x, y, z). Let the point P(X, Y) on the wing plan be enclosed in a small surface 
symmetrical about z=0, the upper half of which will be called S,, and which inter- 
sects the plane z=0 along a closed curve C,. Define a closed surface S$ consisting 
of Sp, the part of the plane z=0 outside C, and a large hemisphere, centre P (X, Y), 
in the upper half space z>0. If dS denotes an element of the surface area of S, 
then Green’s theorem applied to the functions u and|r-R{|~-' harmonic and 
regular in S, gives 
l 


dnJ 
= denoting differentiation along the normal to dS. (It is necessary to exclude 
the point P(X, Y) from the domain of integration by means of Sp as |r—R~' 
is not regular there.) The contribution to the surface integral from the large hemi- 
sphere tends to zero as the radius tends to infinity. By use of the relation 
ou ow 


Haat = and integration by parts, the surface integral is transformed into 


l Ou 


a j w (x, y) (x—X) dxdy f w (x, y) sin 6 dl 


[(x-X¥+(y- YF}? 0. (11) 


The second integral is taken over the part of the wing plan outside C,; in the last 
integral 6 is the angle between the x-axis and the tangent to Cp, d/ is an element of 
arc of C), and the integration is in the positive direction round P (X, Y). 


The second integral may be evaluated on Y=0 by choosing as the curve C, 
segments of the lines x - X +| y | tanG= +e; x- X—| y | cot B= - 2a ccosec 28, where 
a is a small number (Fig. 2). The surface S, consists of sections. of right circular 
cylinders. If e->0 and then a0 the integral over S, tends to —2z u(X, 0), the 
integral over the wing to Ludwieg’s value, and the line integral omitted by Ludwieg to 


2w (X, 0) cos log 
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2é 
= 
2a 
Fig. 2. 
The operation P,. 
Hence 
(X, 0)= +2w (X, O)cos Blog (42) 


Here P, denotes the operation of evaluating the integral by 


Ludwieg’s method, lim lim. 


Alternatively the curve C, may consist of segments of x-X=+e¢, y=+a 
(Fig. 3). In the limit the line integral tends to zero, whence 


w (x, y) dxdy 
(X, 0)=—Px || ° (12a) 


where P, denotes the operation lim lim. 
a0 


For the sake of comparison, consider the case when (X, Y, Z) is not on the 
wing, Z>0. Then P may be surrounded by a‘ small sphere of radius e. When 
e—0 Green’s theorem gives in the limit 


w (x, y)(x—X) dxdy 
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Fig. 3. 
The operation P,. 


Comparison with (12a) shows that when the definition of the surface integral makes 
the line integral + 0, the convergence is non-uniform; 


limu(X, Y, Z)=u(X, Y, 0), .. 
Z-0 
but 


lim w (x, y) (x—X) dxdy {| w (x, y) (x—X) dxdy (14) 


3. KINKED DOUBLY-CYLINDRICAL WINGS 
Ludwieg has considered the case 
w (x, yy=w(x+|y|tanf8) for -l<x+|y|tanB<l 
=0 outside 


where £ is the angle of sweepback. Here the operation P; is appropriate for points 
on the centre line: 
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where the integral is Cauchy’s principal value and 


This relation is due to Neumark, IV.17, and Kiichemann, eqn. 15. (In Kiichemann’s 
report the integral has the wrong sign.) Kiichemann has pointed out that the 
equation may be used to obtain the function w (x) from a given pressure distribution 
u (x), and has indicated a numerical method of calculating w. However, the equation, 
considered as an integral equation for w has an explicit solution (see Appendix) 


an 


= 1 (1 ~ (7). 
w (x)=sec sin cos ral u (x)+ © (1 


(16) 


where = is defined by 
cotz7z=f(8)>0, 0<z<j 


and the condition (see equation (6)) 


1 
w (x)dx=0 
1 
has been used to ensure uniqueness. 


The limiting case «=4 corresponds to a wing with no sweepback, and «=0 
would correspond to the unrealisable case of a sweepback of 90 degrees. When 
x=4, the integral equation (15) and its solution (16) reduce to the well-known 
equations of two-dimensional theory. 


As an example of the use of (16) suppose that 


cos 8} 1+ 1~x | 


x 1-3 
1+x | 


17) 
( 


u(x)= - log 


which gives the pressure distribution at infinity on a wing whose equation is 

C(x, y)=8 y | tan BY J. 
This is the example discussed by Kiichemann in his paper. The integrals arising 


from the substitution of the value of u (x) in the solution may be obtained by contour 
integration, whence 
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nye 
=w (x)=28V ox — 28V, (1 (1 +.x)* cos (22 22? -2x2+x- x’). 


Near x= +1 


~~ — 276 (1— cos m2 . 2? <0. 


— (1 +x)"* cos m.(1—<) <0. 
at 

But >0 at x= 1 
Ox 


is a Necessary geometrical condition. Thus the wing contour obtained in this way 
crosses itself and has no physical reality. (The wing contour is shown crossing 
itself in Kiichemann’s Fig. 16.) ' 


Kiichemann suggests that the design of wings should proceed by successive 
approximation. The wing being supposed to have infinite aspect ratio, the pressure 
distribution at infinity is easily found from an assumed wing contour at infinity. 
The integral equation (15) is then used to find the equation of a wing surface of the 
form z=( (x, y)=((x+| y| tan 8) which will have the same pressure distribution 
along its centre line y=z=0. If a wing is now designed having the given contour 
at infinity, the calculated contour at the centre section, and a smooth transition 
between them, it may be expected that the pressure distribution will be the same 
over every section. The foregoing example shows that this procedure may lead to 
meaningless results when u (x) becomes logarithmically infinite at the edges of the 
wing. 


5 
+ cos l>p>0, >q>0. 


It will be shown that p and q can be chosen in such a way that at infinity and at 
the centre section 


Suppose now that u (x)= 


nye 

==> 0, x= 
ox 
at 
0, x= +1, 
OX 


w (x) <0, 
w (x)>0, x= +1. 


The value at infinity may be deduced from the solution of the integral equation 
for «=4: 


w (x)= 
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| 
If p >4, q > 4, the asymptotic expressions for w (x) near the wing edges are easily 
deduced from 
dx=ny cot ma, O<a<1 
——dx=n 
cot 7a, a< 


Near x= - 1 


w (x) ~- 2-4 = cot x(p—4). (1+) 


Near x= +1 


cot x (g-4). (1-2) a> 0. 


w (x) ~ =y 


The value on the centre line is 


w (x)= sin 7x COS 72 [a +x)? 


tan 7x 


4 


where it is now supposed that q > l—« >}. 


Asymptotic expressions are: 


w (x) ~ sin 7a cos ma. 2-¢[ (1+ x)-? tan cot (p—2).(1+x)?] <0 


near x= 1 if >1, ie. < p<2z, 


and 


near x= +1 if 
The inequalities are satisfied if 


The last equation requires that 2x > 4, 8 < 66.5 degrees. With these restrictions 
p and g must both be less than 1; 0¢/0x will then be integrable over (—1, 1). 


110 


I 
| 

— 

a 


FLOW ROUND SWEPT-BACK WINGS 


This example shows that under suitable conditions the behaviour of the profile 
near the edges is dominated by the behaviour of the pressure near the edges. Too 
much weight should not, therefore, be attached to the occurrence of profiles twisted 
near an edge as it is known from two-dimensional theory that the linear approxi- 
mation is unsatisfactory there. In the above examples the pressure at the stagnation 
point will be finite and large, and not infinite as predicted by linear theory. 


4. EFFECT OF A DISCONTINUITY IN 0¢/¢y ALONG THE CENTRE LINE 


The wing plans considered by Neumark and Kiichemann have discontinuities 
in &X/y across the centre line y=z=0. The effect of rounding off the discontinuity 
will be considered in this section. For the sake of simplicity each half of the wing 
is taken as cylindrical and of constant chord, its equation being 


z={(x, y Dat x+|y | tan B41. 


where g (1)=0, 2’ (1)=0, (0)= tan 8, g(x)=0 for x > 1. is a small parameter. 
X{/ay is now equal to zero along the centre line, except when /=0 (Fig. 4). 
The difference between the velocities u of ¢ (x, y, J) and ¢ (x, y, 0) at a point on 


the centre line, taken to be the origin, is given by Au, where 


xdxdy 


) 


plus the limit of a line integral depending on the method of evaluation of 
the surface integral. 


The integrand of the double integral vanishes except in | y|</. Also, it is 
easy to see that the contribution from the region | x | > a is O (I). 


Now by Taylor’s theorem 


| 
w[xtly tan —~w[x+]|y| tan 6] 


ly |\[ ow ew y | 
1. of iy | tan ¢ ( Ay} +0(x) |, 0<6<1, 


whence the integral over | x | <a, | y| <Jis 


4 [dy dx G+yy? ox ytan8+4l.¢ ( ) +O (I*) 


ay x71. g(| y/l)) 


o (x? + 
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x+|y| tan = CONSTANT 


X = 


<A 


Fig. 4. Rounding off the nose. 
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In this last expression the integration with respect to x can be carried out 
explicitly. The corresponding line integral vanishes. Hence 


1 
— 27Au=40 dt g(t) | log + log + log (a+ - — = ] 


dx’, 


fl 
log dt g(t) —>0 (20) 


as | tends to zero. 


There is thus a smooth transition from a finite large curvature to a 
discontinuous gradient. 


5. WINGS WITH KINKED ISOBARS 


The computations of Neumark and Kiichemann have shown that in the centre 
section the isobars are in general at right angles to the main stream. It may be 
asked whether it is possible to design a wing in such a way that the isobars are 
always inclined at an angle to the main stream and have a kink in the middle. It 
will be shown that unless the surface has a steep ridge over the centre section this 
is not possible without the introduction of special devices, such as suction slots. 


Let 9 be the angle of sweepback of the isobars, © being a function of x not 
zero everywhere. Along each isobar o¢/%x is constant, so that by differentiation 
along the curve 


dy 
_ 
Ox? dydx dx 


But 2 = - cot 9, so 


d 
dydx tan ) 2 
Hence, by integration in the x-direction, 


*?/0x? has only a finite number of zeros. (A and B denote here the limits of 
integration.) Let A, B be chosen in such a way that 0°9/0x? 0 and 9 (x) +0 
in B<x<A. A and B may clearly be arbitrarily near each other. 
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Put A — B=h and consider the potential function 


o(x+h, y, z)- (x, y, 2=9, Y, 2) 


Then at x=B, 


dy = an dx 


is not zero; by continuity there is an interval including B where 0¢,/y 0. It 
follows that there is a cross-flow component on the centre line. It must now be 
examined whether the existence of such a component is compatible 
with linear theory. 


Take the point B as the origin and consider the component v,=9¢,/%y of the 
velocity. It must have a discontinuity across y=0, since v, (X, Y) is an odd function 
of Y. Now 


w, (x, y) Y) dxdy (23) 


- 


(0, Y)= Py 


where the operator P, implies that the first integration is with respect to x and the 
second with respect to y, and w,=99,/%z. To determine the discontinuity in 
v, (0, Y) it is sufficient to take as the domain of integration the rectangle | y | <e, 
|x| <8, and to suppose that 0< | ¥|<e. The contribution from outside the 
rectangle is clearly a continuous function of Y. Choose 4, « so small that 
in the rectangle 


| w, (x, y)—w, 0, y)| <G(y) |x}, 


where G (y) is a bounded function of y. 


Then 
(0, ¥)(y- ¥)dxdy | | x|G 0) | y-¥ | dxdy 
Go)dy 
and 


| (0, y)(y- Y)dxdy _ 4p | y)dy _ 
[x?+(y- 
(0, y)(y- Y) dy 8 


|y-Y |? 
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Hence 


‘2nv, (0, Y)+2P ©. Goray+ 0, y)|dy|y-¥| 
| « «i 


<2 G (y) dy+ =| iw, y)| dy. 
As Y passes through zero, v,(0, Y) has a discontinuity D,, while the 
discontinuity D, in 


* w, 0, y)dy 
y-F 


is easily seen to be independent of «. Thus 
4e 
| 27D,+2D, |< G (y)dy+ lw, @,y)|dy. . 


The left side of the inequality is independent of «, while the right side can be made 
arbitrarily small by choice of «. The left side is therefore zero, i.e. a discontinuity 
D, in v, implies a discontinuity — =D, in 


* w. y)dy 
P| y = 
and conversely. 


It will now be shown that the latter integral is continuous at Y=0. w, (0, y) 
is an even function of y; the difference between the integrals at Y and — Y is 


aw, (0, y)dy 
aly | 


It may be reasonably assumed that w, (0, Y) satisfies a condition 


iw, 0, y)—w, 0, Y)| <M] y-Y|\ y, ¥>0, O<A<1. 


Then 
M |y-Y |*dy 
+2|, ©, log S417 
but du=0 (1) as Y—>0 
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and w, (0, Y) log mas is clearly O (1) since w, (0, Y) is bounded. 


It follows that the integral 


w, (0, y)dy 


is continuous at Y =0, whence it follows that v, (0, Y) is also continuous. But it 
has been seen that there must be a discontinuity in v, if the isobars cross the centre 
line at an angle other than 7/2. This contradiction shows that the isobars must 
cross the centre line at right angles. 


6. UNIQUENESS 


It has just been seen that to a prescribed wing plan A and a given distribution 
of u (e.g. one with kinked isobars) there does not always correspond a profile ¢ (x, y). 
But it is easy to prove that there is either one profile or none. 


Suppose that two profiles exist; let ((x, y) denote their difference and let ¢ 
denote the difference between the corresponding potentials. 


Then 
dep 
—dx=0, z=0,yconstant. ‘ 
Ox 


@ is continuous across the edge of the wing, and | grad¢|=O(R™), O< nu <l, 
where R denotes the distance from the edge. For the sake of simplicity it is assumed 
that A is finite. 


Let S be the closed surface consisting of a large hemisphere in the upper half- 
plane and the plane z=0; the edges of the wing are excluded by small surfaces. By 
one of Green’s formule 


where the volume integral is taken through the interior of S, and = denotes 


differentiation into S. It is easy to see that the contribution to the surface integral 
from the large hemisphere tends to zero as its radius tends to infinity, since 
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/%r=O (r-*). The contribution from the small surfaces used to exclude the edges 


also tends to zero as the surface area tends to zero. In the limit, therefore, 


\\| | grad @ |*dr= - \| dxdy. (28) 
JeA 
But on A, . =0, so ¢=F (y), and also ra dx =0, whence 
This implies 
and the result follows. 
APPENDIX 


An Integral Equation 


The equation connecting the velocity components u (x) and w (x) on the centre 
line is of the form 


: . Gl) 


tanzx{) g(u)du 


1 
where (x) dx=0. 
1 


This equation may be treated by the methods of the theory of functions of a 
complex variable, as was first done by Carleman."*? 


g (u) du 
u-Z 


1 
G (z) 
J-l 
defines an analytic function one-valued in the plane cut from —1 to +1. 
1 
G (x +0) P| (x) .-1x| <1, 
1 


1 
-| 
1 


At infinity 
‘a 
|G (z)|=O(| z|-?), since | ¢ (u) du=0. 
1 
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It is convenient to consider the function 
H (z)=(z+ 1¥ (z- G (z) 
on the real axis y= +0. 


When | x | > 1, H (z) is clearly real, i.e.  H (z)=0. 


When | x | <1, H (z)= (1 - x)'* (Cos —i sin g (u) du 


(x)) 
of which the imaginary part is 


g (u) du 


1 
(1+xy (1 - x)! [sin | 
1 


cos | 
2 
= — 7 cos (1 +x} (1 f (x) 


from the integral equation. The imaginary part of H (z) is therefore given every- 
where on the real axis, and 


| H(z)|=O(|z|~) 
at infinity. 
H (z) is then given everywhere in the upper half-plane by Poisson’s integral 


1 
H | JS H (u) du 
J4 


In particular, on the real axis 


7 1 


(L+uy (1 u)' du 


= — COS | 
But G (z)=(1 - (1+z)-* A (2) 
= —(1—x} (14+. x)-? (Cos 72 +i sin 72) (ZH (z)+i SH (z)) 
of which the imaginary part is = g (x), i.e. 


tan 7 


1 
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Uniqueness 
From the construction it is clear that the formula gives one solution, To show 
that it is the only one the following assumption is sufficient : 
The integrable solution g(x) has a differential coefficient in the open interval 
satisfying 
/ 1 
‘(x)=0 x*—>l. 
(a 
This is satisfied for all functions f (x) occurring in practice. It follows immediately 
that 


g (x)= +5): r—>l. 


Lemma 
log (1 z*) 


It is sufficient to consider the neighbourhood of z= 1. 


G(2)=0| ] near z= +1 


1 
(1-z¥ | (l-u+u zy g (u) du 
1 


u-Z 
= (1 —u)’ g (u) du 
1 u-Z 
Write g, (x)=(1—x¥ g(x) 
=0 |x|>1. J 


g,(x) has by hypothesis a continuous differential coefficient near x=1. Let 
Z=X +iY be any point in the neighbourhood of unity. 


du 


1 
1 u-Z 
du 
+800) 


1 BE | 


1 
(Z)| <| 
The integrand is bounded and | g, (X)| < M, so 
(1- G(Z)=O { log (1-—Z) } 
which proves the lemma. 
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Suppose now that g(x) is the difference between two solutions of the integra! 
equation, so that 


tanza g(u)du 
| 
1 


and g(x) is integrable and O ay) 


Consider the function 
H (x)=(z+ 1) (z—-1)'* G (2). 


The imaginary part of H(z) vanishes everywhere on the real axis. Define 
H (z)=H (z) in the lower half-plane. By Schwarz’s reflection principle, H (z) is 
an analytic function defined in the whole plane, whose only singularities are poles 
or isolated essential singularities at +1. 


But at Z=1 
H (z)=O { log(1-z)} ,0<z2<l, 


so the singularity at z= 1 is a pole of order unity at most and 


a 
H (z)= ——, +a regular function. 


Here a is real since H(z) is real on the real axis. Now as z tends to 1+0 
along the real axis, 


1 1 


1 
| 


and is bounded, i.e. | (z— 1)» H (z)| is bounded along this path. Hence a is zero; 
H (z) is regular at + 1; similarly it is regular at — 1, and therefore in the whole plane. 
Also at infinity H (z)=O(|z|~'). By Liouville’: theorem H (z) is therefore zero 
in the whole plane. Hence g(x)=0 on the real axis, and the solution of the integral 
equation is unique. 


1 
If the condition | g (u) du=0 is dropped, the argument shows that 
1 


tan g (u) du 
1 


(x)= 


has the unique eigenfunction (1+x)-*(1-xy-'. As « goes from 0 to 1, (tan 72)/z 
takes all values from —°© to +0. Thus the integral equation has a continuous 
spectrum consisting of all points on the real axis except 0, for which value the 
equation becomes trivial. 
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Addendum 


For reference some functions occurring in the calculation are given here: 


sin 7x (1+uy (1 —u)'~ u"du 
q u-Xx 


1 


NOTATION 


dS 

Sp 

u, Vv, W 
Vi, W, 
x,y,z 


a 


B 
5, 
(x, y)=eZ (x, y) 
P= 


RAH (2), F H (2) 
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— cos zz (1+ x) (1 — 
=x- 142: n=0 
=x*—x+22 (x - 1)+22? n=1 
— x? +22 (x? - x)4+227 (x- 1) 

3 3 2 


a small parameter. 
wing plan. 
curve of intersection of A and Sp». 
element of arc of Cy). 
in Section 4, a small parameter. 
a fixed number. 
operations of determining principal values of integrals, see 
Figs. 2 and 3. 
numbers satisfying equation (18). 
position vectors. 
element of surface area. 
surface defined in Section 2.1. 
velocity of the main stream. 
components of velocity due to the presence of the wing. 
components of velocity derived from potential ¢,. 
components of position vector r. 
components of position vector R. 
1+sin 8 
defined by = cot zz=log 
angle of sweepback, see Fig. 1. 
arbitrarily small constants. 
equation of wing surface. 
angle of sweepback of isobars. 
density. 
velocity potential. 
defined by equation (22). 
real and imaginary parts of the complex function H (z). 
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